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ABSTRACT 


In  1962,  the  United  States  eondueted  its  final  atmospherie  nuelear  test.  Sinee  1962,  the 
Ameriean  national  laboratories  have  attempted  to  simulate  the  results  observed  in  exo- 
atmospherie  testing  in  order  to  understand  and  explain  how  high-energy  eleetrons  beeame 
trapped  further  aeross  the  Earth’s  magnetie  field  than  expeeted. 

In  this  thesis,  my  researeh  will  use  a  eomputer  modeling  program  designed  by  the 
late  Dr.  Dennis  Hewett  for  LLNL  in  1972  to  simulate  the  debris  from  a  High  Altitude 
Nuelear  Explosion  (HANE).  The  objeetive  of  this  researeh  is  to  examine  two  physieal 
phenomena,  eollisionless  shoeks  driven  by  multiple  ion  speeies  and  eharge  exehange,  to 
determine  their  independent  relevanee  to  the  final  spatial  disposition  of  fission  fragments 
from  a  HANE. 

This  researeh  used  the  ZMR  eode,  a  one-dimensional,  partiele-in-eell  plasma  eode 
to  simulate  the  movement  of  the  debris  ions  produeed  by  the  HANE.  The  debris  ions  are 
assumed  to  be  the  souree  of  the  measured  high-energy  eleetrons  due  to  subsequent  beta- 
deeay.  These  high-energy  eleetrons  ean  damage  the  satellite  network  vital  to  the 
Department  of  Defense  and  seeurity  of  the  nation. 
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I.  INTRODUCTION 


A.  INTRODUCTION 

STARFISH  PRIME  (hereafter  referred  to  as  STARFISH)  was  a  high-altitude 
nuelear  explosion  (HANE)  condueted  by  the  United  States  (U.S.)  on  9  July  1962.  While 
the  explosion  oeeurred  at  400  kilometers  (km)  above  the  Earth’s  surfaee,  satellites 
measured  high-energy  fission-spectrum  electrons  at  locations  (elevations)  and 
corresponding  magnetic  field  lines  well  above  what  analysts  have  been  able  to  explain. 
Bomb  debris  from  STARFISH  was  composed  of  fission  products,  other  ionized 
materials,  and  electrons,  some  born  from  beta  decay  of  fission  fragments.  Charged 
particles  moving  in  the  presence  of  the  existing  magnetic  field  of  the  Earth  will  gyrate 
around  the  field  lines.  The  binding  of  these  charged  particles  to  the  field  lines  limits  how 
far  charged  particles  will  travel  in  the  direction  perpendicular  to  the  field  lines.  That  is 
part  of  the  mystery  observed  with  the  STARFISH  test;  the  question  of  how  ions  and 
electrons  were  transported  to  higher  field  lines  has  not  been  fully  explained  by  analyses 
using  simulations  relying  on  fluid-based  principles  suggesting  that  kinetic  effects  and 
other  physical  processes  must  be  incorporated  to  develop  a  full  description  of  the  test. 

The  measured  energetic  electrons  from  STARFISH  are  more  than  just  a  scientific 
curiosity.  Electrons  and  ions  from  a  HANE  can  be  expected  to  collide  with  satellites  and 
the  resulting  accumulation  of  damage  can  result  in  a  corresponding  increase  in  failure 
rates  (i.e.,  a  reduction  in  the  mean-time-to-failure).  Following  STARFISH,  three  satellites 
quickly  failed  and  another  two  failed  well  before  their  predicted  lifetimes.  The  dramatic 
increase  in  the  number  and  importance  of  military  as  well  as  civilian  satellites  from  1962 
to  the  present  underscores  how  critical  this  effect  can  be  to  the  Department  of  Defense 
(DOD)  and  to  the  Nation. 

In  addition  to  the  anomalous  distribution  of  debris,  the  energy  of  STARFISH 
“dragged”  the  Earth’s  magnetic  field  lines  outward  and  created  a  diamagnetic  cavity  or 
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“bubble.”!  In  such  a  bubble,  the  magnetic  field  strength  is  well  below  that  of  the  ambient 
field,  and  in  this  region  charged  particles  will  experience  a  reduced  level  of  turning  and 
binding  to  the  magnetic  field  lines.  Therefore,  charged  particles  in  a  diamagnetic  cavity 
will  travel  in  a  generally  straight  line  until  they  depart  the  bubble.  However,  the 
diamagnetic  cavity  observed  in  STARFISH  was  not  large  enough  to  fully  explain  the 
spatial  distribution  of  energetic  electrons  measured.^  Additionally,  attenuation  of  the 
cosmic  background  radio  waves  gives  insight  to  ionization  of  the  atmosphere.  There  were 
differences  in  the  attenuation  data  when  compared  from  collection  stations  in  Alaska  to 
collection  stations  in  Canada  and  Norway  leading  researches  to  conclude  “that  upward- 
traveling  debris  from  the  explosion  managed  to  work  its  way  through  the  Earth’s 
magnetic  field  vertically  overhead  of  Johnston  Island.”^ 

Fission  of  an  actinide  generally  produces  two  ionized  fission  fragments.  These 
fission  fragments  are  highly  ionized;  the  process  of  fission  injects  both  positive  ions  and 
negative  electrons  into  the  background  plasma.  Additionally,  these  fission  fragments  are 
normally  neutron  rich  and  decrease  their  internal  (nuclear)  energy  by  beta  decay. 
Generally  a  fission  fragment  undergoes  several  beta  decays  before  becoming  a  stable 
isotope.^  These  beta  decays  are  a  source  of  energetic  electrons.  It  is  unknown  where  the 
energetic  electrons  measured  after  STARFISH  came  from;  however,  one  possible 
explanation  (based  on  the  measured  energy  spectrum  and  the  fact  that  they  were  not 
present  before  the  detonation)  is  they  were  from  the  beta  decays  of  fission  fragments  that 
had  been  somehow  transported  to  much  higher  altitudes  than  their  point  of  generation. 

This  work  considers  two  physical  processes  as  possible  explanations  for  the 
presence  of  the  large  number  of  energetic  electrons  at  field  lines  above  where  they  were 
expected.  The  two  processes  explored  here  are  a  multispecies  collisionless  shock  and 

!  Palmer  Dyal,  “Particle  and  Field  Measurements  of  the  STARFISH  Diamagnetic  Cavity,”  Journal  of 
Geophysical  Research  111  (2006):  1. 

2  Ibid.,  13. 

2  R.P  Balser,  R.B.  Dyce  and  H.  Leinback.  “High-Latitude  Ionization  Associated  with  the  July  9 
Explosion,”  Journal  of  Geophysical  Research  68  (1963):  741. 

^  K.Krane.  Introductory  Nuclear  Physics  (New  York:  John  Wiley  &  Sons,  1988),  485. 
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collisional  charge  exchange.  A  collisionless  shoek  oceurs  when  the  expanding  plasma 
generates  an  electric  field  of  large  enough  intensity  to  rapidly  accelerate  (shoek)  the 
background  plasma.  When  the  expanding  plasma  consists  of  more  than  a  single  ionic 
constituent,  it  is  ealled  multispecies.  Is  it  possible  a  debris  species  with  a  small  gyro- 
radius  would  form  a  collisionless  shock  that  would  subsequently  accelerate  other  speeies 
with  larger  gyro-radii?  Charge  exchange  is  the  proeess  where  one  (or  more)  eleetron(s)  is 
exehanged  between  atoms.  The  change  in  charge  state  affects  how  the  ion’s  trajeetory 
bends  in  the  magnetic  field.  Neutral  atoms’  trajectories  do  not  bend  in  the  magnetic  field 
and,  they  travel  in  a  straight  line  until  a  eollision  diverts  their  path.  Could  this  process 
explain  the  electron  population  further  from  the  blast  point  and  aeross  many  magnetic 
field  lines?  Those  two  questions  are  addressed  in  this  work. 

These  two  processes  were  modeled  using  a  cylindrical  geometry  eode  called 
ZMR.  ZMR  was  initially  developed  (and  subsequently  updated)  at  Lawrenee  Livermore 
National  Laboratory  (LLNL)  30  years  ago  and  has  been  used  to  simulate  various  plasma 
phenomenologies  principally  in  cases  where  the  motion  of  the  ions  was  important  but  the 
detailed  eleetron  dynamies  do  not  influence  the  phenomena  of  interest.^  This  situation 
arises  when  the  Debye  length  is  small  eompared  to  the  ion  gyro-radius,  so  the  quasi¬ 
neutral  assumption  is  valid  and  the  electron  inertial  length  is  much  smaller  than  the  ion 
inertial  length.^  This  dissertation  researeh  expanded  upon  the  multispeeies  investigation 
began  by  Hewett  and  then  added  a  eapability  to  examine  charge  exchange.  This  work  is 
eomputational  in  nature  and  relies  on  simulation  to  provide  insight  into  physieal 
proeesses.  Understanding  these  two  processes  in  a  basie  one-dimensional  (1-D)  geometry 
will  help  in  the  analysis  of  their  importance  and  may  lead  to  their  incorporation  into  more 
robust  future  simulations.  That  is  the  goal  of  this  work.  This  research  is  not  intended  to 
provide  a  full  reproduetion  the  STARFISH  data  set.  However,  this  work  is  a  preliminary 


^  D.W.  Hewett.  “A  Global  Method  of  Solving  the  Electron-Field  Equations  in  a  Zero-Inertia-Electron- 
Hybrid  Plasma  Simulation  Code,”  J.  Comp.  Phys.  38  (1980)  378-395. 

®  G.F.  Simonson,  D.W.  Hewett,  DE.  Shumaker,  D.V.  Anderson,  J.A.  Byers  and  G.DiPeso.  “Numerical 
Simulation  of  Interactions  of  Nuclear  Weapons  with  the  Geomagnetic  Field,”  UCRL-JC-109841 :  Sensor 
Technical  Conference  (1992),  5. 
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investigation  of  relevant  physieal  hypotheses  that  eould  aeeount  for  observed  anomalies. 
From  this  work,  it  will  be  possible  to  integrate  and  test  the  developed  eomputational 
methods  into  larger  and  more  robust  HANE  simulation  eodes. 

LLNL,  other  national  laboratories  and  the  DOD  have  been  studying  STARFISH 
sinee  1962  and  attempting  to  explain  the  results.  The  limited  number  of  high  altitude  tests 
and  the  limited  data  from  eaeh  test  hinder  the  eommunity’s  overall  understanding  of  the 
eomplex  problem.  Simulation  models  will  inerease  in  eomplexity  and  fidelity  as 
eomputing  power  inereases.  Original  eodes  used  magneto-hydrodynamies  (MHD)-based 
modeling  where  the  plasma  was  treated  as  a  eharged  fluid  where  the  rules  of  fluid 
dynamies  are  eombined  with  the  rules  of  eleetrieity  and  magnetism.  Subsequent  inereases 
in  eomputing  power  allowed  the  use  of  the  Partiele-In-Cell  (PIC)  method.  PIC  eodes  are 
used  when  the  assumption  of  a  Maxwellian  veloeity  distribution  is  inadequate  to  deseribe 
the  plasma.^  At  the  altitude  of  STARFISH,  the  mean  free  path  of  the  debris  ions  in  the 
thin  atmosphere  is  too  large  to  use  a  MHD  eode.^  ZMR  is  a  PIC  eode  with  an  implicit 
field  solver^  and  it  provides  better  resolution  for  the  spatial  motion  of  the  ions  than  an 
MHD  code. 

B,  PREVIOUS  WORK  IN  THIS  FIELD 

Since  1962  there  have  been  many  varied  efforts  to  explain  the  STARFISH  test. 
The  Journal  of  Geophysical  Research  dedicated  its  whole  68th  volume  to  the  analysis 
and  presentation  of  the  data.  The  editors  stated  in  the  introduction  that,  “a  central  purpose 
of  this  introduction  is  to  call  attention  to  the  fact  that  there  still  exists  a  lack  of  agreement 
on  the  spatial  extent  of  the  STARFISH  radiation  belt  and  on  the  total  number  of  trapped 
particles  therein.”!*'  The  journal  devoted  that  issue  to  the  efforts  of  various  authors  and 


^  E.C.  Whipple  Jr.  “A  kinetic  approach  to  magnetospheric  modeling,”  Quantitative  Modeling  of 
Magnetospheric  Processes  (American  Geophysical  Union:  Washington,  DC,  1979),  463. 

^  T.G.  Cowling,  Magnetohydrodynamics  (Bristol  UK:  Adam  Hilger,  1976):  120. 

D.W.  Hewett,  Documentation  of  the  EM  2.5D  Quasi-neutral  Hybrid  Code  ZEMER,  10. 

!**  W.L.  Brown,  W.N.  Hess  and  J.A.  Van  Allen,  “Collected  Papers  of  the  Artificial  Radiation  Belt  from 
the  July  9,1962,  Nuclear  Detonation,”  Journal  of  Geophysical  Research  68  (1  February  1963):  605. 
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institutes  to  capture,  analyze  and  present  the  data  on  the  introduced  radiation  measured 
following  STARFISH. 

The  Trapped  Radiation  Handbook  was  originally  compiled  in  1971  for  the 
Defense  Nuclear  Agency  and  was  revised  several  times  in  the  intervening  years. 
Originally  designed  as  a  single  source  for  satellite  proponents  to  understand  the  trapped 
radiation  environment,  it  evolved  into  a  repository  of  derivations  and  facts  related  to 
charged  particles  in  the  Earth’s  magnetic  field. The  handbook  provided  facts  and  data 
about  trapped  radiation  but  did  not  attempt  to  establish  from  first  principles  how  the 
artificial  radiation  became  trapped  in  the  radiation  belts.  While  the  experts  have  laid  out 
possible  mechanisms,  none  were  conclusively  examined. 

In  1978,  Kilb  and  Glenn  employed  a  collisionless  magneto-hydrodynamic 
(CMHD)  simulation  that  included  charge  exchange  for  the  atmospheric  ions.^^  For 
STARFISH  and  the  relatively  low  density  of  the  atmosphere,  the  authors  claimed  the 
fluid  equations  of  MHD  did  not  produce  the  observed  results. However,  the  authors 
use  a  CMHD  model  below  the  burst  to  specifically  look  at  the  formation  of  a  low  altitude 
ion  patch.  Additionally,  their  modeling  employed  charge  exchange  but  only  for  the 
atmospheric  constituents  of  oxygen  and  nitrogen;  the  debris  ions  were  not  modeled  in  the 
CMHD  simulation.  16  In  the  simulation,  the  atmospheric  ions  were  neutralized  and  then 
moved  in  a  straight  line  across  magnetic  field  lines.  At  a  later  time,  the  neutral  air  atoms 
were  re-ionized  to  form  the  ion  patch,  i^  The  modeling  results  appeared  to  conserve 
energy.  1^  The  model  focused  on  the  total  amount  and  spectra  of  the  ultraviolet  radiation 

1 1  Trapped  Radiation  Handbook. 

12  Ibid 

12  R.W.  Kilb,  and  D.E.  Glenn.  CMHD  Simulations  of  very  high  altitude  nuclear  bursts  from  0  to  I 
second  (U)  (Santa  Barbara:  Mission  Research  Corporation,  1978),  17 

l^Ibid.,  17. 

15  Ibid.,  17. 

16  Ibid.,  38. 

12  Ibid.,  37. 

18  Ibid.,  61-62. 
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generated^^  and  not  the  spatial  distribution  of  the  debris.  The  ehoice  of  a  CMHD  model 
prevented  the  investigation  of  eollisionless  shock  as  this  is  a  kinetic  ion  driven  process. 
Additionally,  Kilb  and  Glenn  focused  on  charge  exchange  for  constituents  of  the 
atmosphere  only  and  did  not  examine  the  debris  ions.  As  will  be  discussed  in  Chapter 
II. C,  the  debris  ions  likely  undergo  several  beta  decays  and  eject  high-energy  electrons, 
which  contribute  to  the  electron  population  measured  following  STARFISH. 

LLNL  published  simulations  using  a  3D  hybrid  particle  code  originally  designed 
for  theta  pinch  plasma.^o  Shumaker  used  a  3D  hybrid  code  called  QN3D  for  3D 
simulations  of  a  single  debris  species  expanding  into  a  background  plasma.  The  objective 
of  that  research  involved  investigation  of  the  formation  of  the  diamagnetic  cavity. 
Additionally,  by  collapsing  a  spatial  dimension,  the  code  produced  2D  results  for  the 
study  of  longer  time  kinetic  effects.^^  However,  the  code  required  long  run  times  with 
only  a  single  debris  species^^  and  is  no  longer  in  use  at  LLNL. 

The  review  of  the  literature  did  not  identify  any  work  that  examined  either 
collisionless  shock  driven  by  multiple  ion  species  or  charge  exchange  as  processes  to 
move  debris  ions  across  magnetic  field  lines.  These  debris  ions  would  be  the  source  of 
high-energy  electrons  by  beta  decay.  Additionally,  these  processes  may  be  important  to 
other  simulation  work  at  LLNL  and  elsewhere. 

C.  RESEARCH  OBJECTIVES  AND  METHODOLOGY 

This  research  in  this  dissertation  focused  on  two  physical  processes — collisionless 
shock  driven  by  multiple  ion  species  and  charge  exchange.  These  processes  were 
examined  in  relation  to  the  debris  from  a  HANE  expanding  into  the  Earth’s  ambient 
magnetic  field  and  upper  atmosphere.  They  were  examined  to  judge  their  relative  impact 

19  Ibid.,  58. 

E.  Shumaker,  D.  V.  Anderson  and  G.F.  Simonson.  “3D  and  2D  Simulations  of  an  Expanding  Plasma 
Using  a  Darwin  Hybrid  Particle  Code.”  UCRL-ID-109852:  LLNL:  1. 

21  Ibid.  10. 

22  Ibid.  1. 

23  Ibid.  12. 
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on  the  spatial  distribution  of  the  debris  following  a  HANE.  This  work  foeused  on  the 
debris  from  a  HANE  sinee  the  debris  is  a  souree  of  high-energy  upon  eleetrons  following 
beta  deeay.  High-energy  electrons  were  measured  far  across  the  Earth’s  magnetic  field 
lines  following  a  high  altitude  nuclear  test.  These  high-energy  electrons  can  impact  on 
satellites  and  damage  their  electronics,  thereby  decreasing  the  mean-time-to-failure.  The 
U.S.  Government  relies  heavily  on  satellites,  especially  during  a  crisis.  Our  continued 
reliance  on  satellites  makes  us  vulnerable  to  the  effects  of  a  high  altitude  nuclear 
explosion.  Any  improvement  in  our  ability  to  model  and  forecast  the  spatial  distribution 
of  high-energy  electrons  following  a  HANE  will  increase  our  ability  to  protect  and 
operate  our  vital  satellite  fleet. 

The  method  of  examination  consisted  of  simulation  using  a  one  dimensional 
hybrid  plasma  code  described  in  detail  in  Chapter  VI.  Simulations  were  performed  in 
simple  radial  geometry  (on  a  faster  running  code  than  a  full  three  dimensional  code)  to 
develop  an  understanding  of  the  importance  of  the  two  processes.  If  either  or  both 
processes  are  judged  to  contribute  to  the  population  of  high-energy  electrons  beyond 
what  has  been  previously  predicted,  then  those  processes  should  be  included  in  a  higher 
fidelity  simulation. 

The  flow  of  the  dissertation  begins  with  an  overview  of  the  primary  physics  topics 
that  are  required  for  this  research.  The  physical  details  of  the  STAREISH  PRIME  nuclear 
test  are  presented.  Then  the  mechanics  of  collisionless  shock  formation  are  derived,  and 
an  overview  of  collisional  charge  exchange  is  provided.  The  final  piece  of  background  is 
an  overview  of  the  hybrid  plasma  code  ZMR  as  well  as  a  discussion  of  the  modifications 
made  to  the  code  for  this  work.  Einally,  the  results  and  conclusions  from  the  research  are 
provided  and  their  implications  discussed. 


7 


THIS  PAGE  INTENTIONALLY  LEET  BLANK 


8 


II.  BACKGROUND  PHYSICS 


This  chapter  explains  in  more  detail  some  of  the  baekground  physies  required  for 
this  work.  Readers  with  a  thorough  understanding  of  the  Earth’s  magnetie  field,  motion 
of  a  eharged  partiele  in  a  eonstant  magnetie  field  and  nuelear  fission  ean  skip  ahead  to 
Chapter  III.  This  ehapter  provides  a  brief  overview  of  those  three  topies  as  baekground  to 
assist  in  understanding  the  remainder  of  the  work. 

A,  L-SHELL  MODEL  OF  EARTH’S  MAGENTIC  FIELD 

The  Earth’s  magnetic  field  resembles  a  simple  dipole  field  and  is  frequently 
modeled  as  sueh.  In  a  dipole  magnetie  field,  there  are  two  opposite  poles  with  magnetie 
field  lines  that  follow  a  eurved  path  between  them.  Eigure  1  is  a  sketeh  of  the  standard 
dipole  representation  of  the  Earth’s  magnetie  field.  A  slight  eomplieation  in  the 
representation  is  that  the  southern  magnetie  pole  is  in  the  Earth’s  northern  hemisphere.^'^ 

Transferring  this  dipole  model  of  the  Earth’s  magnetie  field  into  a  mathematieally 
usable  form  for  modeling  and  simulation  requires  its  representation  in  an  appropriate 
eoordinate  system.  Eor  modeling  the  Earth’s  field,  either  a  Cartesian  or  a  spherieal 
eoordinate  system  eould  work;  however,  eaeh  of  these  representations  has  its  own 
ehallenges.  The  Cartesian  coordinates  system  is  mathematically  easier  to  work  with  for 
vector  equations.  However,  the  Cartesian  system  is  eumbersome  and  intuitively  difficult 
when  used  to  deseribe  loeations  in  a  spherieal  frame  of  reference.  The  spherieal 
eoordinate  system  has  a  more  appropriate  geometry  for  the  Earth’s  magnetie  field  but  is 
eumbersome  to  work  with  mathematieally.  Also,  the  spherieal  eoordinate  system  has  a 
simpler  deseription  for  the  distanee  of  a  point  from  the  eenter  of  the  system,  but  a  more 
diffieult  deseription  for  the  distanee  between  two  points  in  spaee. 

Mellwain  proposed  labeling  a  drift  shell,  whieh  is  a  constant  magnetic  field  line 
rotated  around  the  dipole,  by  two  parameters.  He  used  Bm  as  the  sealar  value  of  the 


J.W.  Chamberlain  On  Motion  of  Charged  Particles  in  the  Earth ’s  Magnetic  Field  (New  York: 
Gordon  and  Breach,  1964),  2. 
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magnetic  field  at  the  mirror  point  and  L,  which  is  generally  the  number  of  Earth  radii 
measured  at  the  magnetic  equator  of  the  dipole. More  specifically 


Re  sin' 0  (1) 

where  r  is  the  radial  distance  from  the  center  of  the  Earth,  Re  is  the  radius  of  the  Earth 
and  0  is  the  angle  measured  down  from  the  pole  of  the  dipole.  While  this  system  was 
originally  designed  to  eliminate  the  need  for  a  numerical  solution  to  a  differential 
equation,  the  system  has  continued  to  be  used  broadly  by  researchers  studying  space. ^6 
The  E-shell  method  is  used  in  much  of  the  literature  describing  STAREISH’s  post 
detonation  data  collection. 


Eigure  1.  Simple  drawing  of  Earth’s  E-Shell  organization  and  general  sketch  of  a 

magnetic  dipole  field. 


Carl  E.  Mcllwain,  “Coordinates  for  Mapping  the  Distribution  of  Magnetically  Trapped  Particles,” 
Journal  of  Geophysical  Research  66  (1961):  3681. 

2®  M.  Schulz,  “Canonical  Coordinates  for  Radiation-Belt  Modeling,”  Geophysical  Monograph,  {1996): 

157. 
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B,  MOTION  OF  A  CHARGED  PARTICLE  IN  A  CONSTANT  MAGNETIC 
FIELD 


The  motion  of  a  charged  particle  subjected  to  electric  and  magnetic  fields  begins 
with  a  description  of  the  Lorentz  force 


F  = 


m- 


dv 

dt 


+  VX  — 


(2) 


where  m  is  the  mass  of  the  particle,  t  is  the  time  element,  q  is  the  charge  state  of  the 
charged  particle,  E  is  the  electric  field,  v  is  the  velocity  of  the  particle,  B  is  the  magnetic 
field  and  c  is  the  speed  of  light  in  a  vacuum.^^ 

Assuming  there  is  no  electric  field  this  becomes 


dv 

dt 


1. 

m 


-  B 
vx  — 
c 


(3) 


Further  assuming  that  the  local  magnetic  field  exists  only  in  the  z  direction  of  the 
Cartesian  coordinate  system,  this  becomes 


dt  me  ^  ‘ 


(4) 


dt  me 


(5) 


dt 


=  0 


(6) 


Thus,  the  velocity  in  the  z  direction  is  constant.  Since  there  is  no  electric  field  and  the 
magnetic  field  does  not  change  the  energy  of  the  charged  particle,  the  initial  z  component 
of  the  velocity  is  unchanged.  Taking  another  time  derivative  of  equations  4  and  5  and 


J.D.  Jackson.  Classical  Electrodynamics  (Hoboken  NJ  :  John  Wiley  &  Sons,  1999),  553. 
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combining  with  the  first  derivatives  of  equations  5  and  4  respeetively,  the  differential 
equations  beeome 


me  dt  \mc  ^ 


-qB^  dv^ 


dr  me  dt 


\  me  J 


These  equations  have  the  standard  solutions  of 


v^(0  =  Cj  sin 


V  me 


t  +  8. 


V  rO  =  C2sin 


V  me 


t  + 


(V) 

(8) 


(9) 

(10) 


C  and  5  are  eonstants  of  the  integration  and  do  not  affeet  the  geometrie  nature  of  the 
motion;  they  do  have  partieular  values  for  a  given  set  of  initial  eonditions  of  the  eharged 
partiele.  Therefore,  the  resulting  solutions  to  the  veetor  eomponent  equations  of  motion 
for  a  eharged  partiele  in  magnetie  field  generate  helieal  motion  along  the  lines  of 
magnetie  force.  It  is  important  to  note  that  the  presenee  of  a  magnetie  field  does  not 
impaet  the  component  of  motion  of  the  eharged  partiele  parallel  to  the  magnetie  field 
(assumed  to  be  in  the  z  direetion  for  this  explanation).  The  term  modifying  the  time 
eomponent  in  equations  9  and  10  is  defined  as  the  eyelotron  frequeney,  tOc. 


me  _ 


(11) 


The  eonstant  of  integration  C  is  the  eomponent  of  the  veloeity  in  the  plane  perpendicular 
to  the  magnetie  field.  This  has  the  symbol  Vj^.  In  this  example,  the  magnetie  field  is  in  the 
z  direetion  so  Vj^  resides  in  the  x-y  plane.  Integrating  equations  9  and  10  again,  the 
equations  for  the  positions  inx  andy  beeome 


12 


and 


(12) 


me 

x  =  x„  +  v,  — cos(cot  + 

°  ^  qB  " 


<5:) 


3'  =  3^o  +  ^i— cos(ft»,?  +  (5,) 

qB 


(13) 


The  term  modifying  the  trigonometry  functions  above  is  defined  to  be  the  Larmor  radius, 

n- 


mcvj_ 

qB  (14) 


The  Larmor  radius  describes  the  perpendicular  radius  of  the  helical  motion  of  the  charged 
particle .28  For  a  constant  initial  particle  velocity  and  magnetic  field,  the  ratio  of  mass  to 
charge  state  determines  the  Larmor  radius.  The  Larmor  radius  helps  describe  how  the 
path  of  a  charged  particle  turns  while  moving  through  a  magnetic  field. 

C.  MOTION  OF  A  CHARGED  PARTICLE  IN  A  DYNAMIC  MAGNETIC 

FIELD 

In  addition  to  the  above-described  helical  motion  of  charged  particles  bound  to 
magnetic  field  lines,  the  phenomenon  of  drift  also  occurs.  There  are  many  physical 
causes  for  such  drifts  and  they  include  the  presence  of  an  electric  field,  non-uniformity 
(i.e.,  curvature  and  gradient)  of  the  magnetic  field  and  the  time  variability  of  the  magnetic 
field.29  Charged  particles  with  opposite  charges  drift  in  opposite  directions.  The  motion 
of  ions  is  explicitly  simulated  in  the  ZMR  code,  so  all  the  relevant  drift  motions  are 
accounted  for. 


28  Francis  Chen.  Introduction  to  Plasma  Physics  and  Controlled  Fusion  (New  York:  Springer,  2006), 

20. 

29  Ibid., 43. 
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Figure  2  shows  the  motion  of  a  eharged  partiele  bound  to  a  magnetie  field  line  on 
the  Earth.  The  partiele  has  a  helieal  motion  as  it  travels  up  and  down  the  field  line.  There 
are  no  non-eonservative  forees,  so  the  particle  does  not  gain  or  lose  energy.  The  portion 
of  the  velocity  that  is  perpendicular  to  the  magnetic  field  is  a  minimum  at  the  magnetic 
equator,  where  the  field  strength  is  lowest,  and  increases  to  100%  at  the  ends  the  helical 
path  where  the  field  strength  is  highest.  The  end  points  are  called  the  mirror  points.  The 
whole  helical  path  of  the  charged  particle  then  drifts  according  to  the  parameters  of  the 
fields  and  any  forces.  Figure  2  shows  the  helical  path  of  the  charged  particle  between  the 
mirror  points  as  well  as  the  drifting  of  the  whole  helical  path. 


Figure  2.  Motion  of  a  charged  particle  in  the  Earth’s  dipole  magnetic  field  (after 

Hess,  1964). 

Absent  any  input  or  outflow  of  energy,  the  charged  particle  will  maintain  this 
helical,  drifting  motion.  There  are  adiabatic  invariants  for  the  motion  of  charged  particles 
around  the  lines  of  magnetic  field  as  well  as  along  the  line  of  magnetic  field  and  the  drift 
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motion  perpendicular  to  the  magnetic  field.^o  The  adiabatic  invariants  are  well 
documented  and  can  be  solved  explicitly.  The  importance  of  the  adiabatic  invariants 
comes  from  the  ability  to  predict  the  motion  of  charged  particles  as  long  as  the  conditions 
remain  constant.  In  such  cases,  the  charged  particles  will  stay  bound  to  the  magnetic 
field.  A  HANE  release  a  large  amount  of  energy  in  a  short  time  and  a  compact  space,  and 
the  subsequent  additional  energy  release  (e.g.,  from  beta  decay  of  fission  fragments)  is 
very  small  in  relation.  Thus,  there  is  no  mechanism  to  input  significant  additional  energy 
into  the  plasma  to  change  the  adiabatic  invariants. 

D.  OVERVIEW  OF  FISSION 

Fission  is  a  nuclear  process  where  a  nucleus  splits  into  two  or  more  smaller  nuclei 
with  an  associated  release  of  energy.  Fission  is  normally  triggered  when  a  nucleus 
(commonly  uranium-235  or  plutonium-239)  absorbs  a  neutron  and  forms  a  compound 
nucleus.  That  compound  nucleus  is  unstable  and  splits  into  two  or  more  fission 
fragments,  although  statistically  the  probability  of  three  fission  fragments  is  small  (1  in 
300).32  The  sum  of  the  masses  of  these  fission  fragments  is  less  than  the  original  mass  of 
the  parent  compound  nucleus  and  the  missing  mass  is  released  as  energy  on  the  order  of 
200  mega-electron-volts  (MeV)  per  fission.33 

The  two  (or  more)  fission  fragments  released  from  each  fission  reaction  do  not 
have  specified  atomic  or  mass  numbers.  Rather,  they  follow  a  bimodal  statistical 
distribution  based  on  mass  numbers  of  the  two  fragments.  Based  on  the  distribution,  the 
products  are  referred  to  as  the  light  group  with  mass  number  between  80-110  and  the 
heavy  group  with  mass  number  between  130  and  155.  The  distribution  of  the  mass  of  the 


30  George  Parks,  Physics  of  Space  Plasmas,  An  Introduction,  2nd  ed.  (Boulder:  Westview  Press, 
2004),  135. 

3 1  Charles  Bridgman,  Introduction  to  the  Physics  of  Nuclear  Weapons  Effects  (Wright  Patterson  AF 
Base:  Defense  Threat  Reduction  Agency,  2001),  77. 

33  H.D.  Bush.  Atomic  and  Nuclear  Physics:  Theoretical  Principles  (Englewood  Cliffs:  Prentice-Hall 
Inc,  1962),  136. 

33  Charles  Bridgman,  Introduction  to  the  Physics  of  Nuclear  Weapons  Effects  (Wright  Patterson  AF 
Base:  Defense  Threat  Reduction  Agency,  2001),  77. 
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fission  fragments  for  fission  of  uranium-235  is  shown  in  Figure  3.  The  distribution  for 
fission  of  plutonium-239  is  similar.  Regardless  if  the  parent  nucleus  is  uranium  or 
plutonium,  both  the  fission  fragments  have  excess  neutron-to-proton  ratios  and  therefore 
frequently  seek  a  more  stable  nuclear  state  by  one  or  more  beta  decays.  Beta  decay 
converts  a  neutron  into  a  proton  and  releases  an  electron  of  variable  energy. 34  The  half- 
life  for  beta  decay  varies  by  isotope  and  can  range  from  milli-seconds  to  years  (with  tens 
of  years  being  the  most  abundant).  Around  4%  of  the  total  energy  released  from  fission 
occurs  from  the  beta  decay  of  the  fission  fragments. 35  This  research  considers 
mechanisms  by  which  fission  fragments  can  be  transported  to  higher  L-shells  where  their 
subsequent  beta  decay  will  introduce  high-energy  electrons  that  may  account  for 
otherwise  unexplained  anomalous  radiation  levels. 


Mass  number  [amu] 


34  H.D.  Bush.  Atomic  and  Nuclear  Physics:  Theoretical  Principles  (Englewood  Cliffs:  Prentice-Hall 
Inc,  1962),  137. 

35  Charles  Bridgman.  Introduction  to  the  Physics  of  Nuclear  Weapons  Effects  (Wright  Patterson  AF 
Base:  Defense  Threat  Reduction  Agency,  2001),  77. 
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Figure  3.  Distribution  of  mass  of  fission  fragments  from  the  fission  of  uranium-235 

by  a  thermal  neutron  (after  Krane,  1988) 

The  energy  speetrum  of  the  eleetrons  released  following  the  fission  of  uranium- 
235  was  experimentally  measured  and  fits  the  distribution 

T  =  3.88exp[-0.575£-0.055£"]  (-^5) 

where  Y  is  the  distribution  [electrons  per  fission  per  MeV]  and  E  is  the  electron  energy 
[MeV];  formula  15  is  valid  for  E  <  7  MeV7®  Energy  dependent  detectors  allow  the 
comparison  of  sampled  electrons  to  the  distribution  of  equation  15  and  a  determination  of 
the  fit  of  the  measured  spectrum  to  the  experimentally  derived  spectrum  of  electron 
energy  from  beta  decay. 

A  nuclear  explosion  consists  of  a  very  large  number  of  atoms  undergoing  fission 
releasing  a  corresponding  amount  of  energy  in  a  very  short  amount  of  time.  The  total 
energy  released  will  raise  the  temperature  at  the  point  of  the  explosion  to  tens  of  millions 
of  degrees  Eahrenheit.^^  The  high  temperatures  achieved  by  a  nuclear  explosive 
detonation  will  cause  a  high  degree  of  ionization  of  the  resulting  debris.  The  ionized 
debris  will  consist  of  both  positive  ions  and  free  electrons.  However,  these  initial  ions  and 
free  electrons  will  be  bound  to  the  Earth’s  magnetic  field  lines  as  they  stream  away  from 
the  burst  point.  Even  the  creation  of  a  diamagnetic  bubble  (as  described  in  the  next 
chapter)  does  not  explain  the  movement  of  the  ionized  fission  fragments  and  their 
electrons  to  elevations  observed  in  STAREISH.  Therefore,  the  initial  production  of 
electrons  is  not  believed  to  be  the  source  of  the  high  E-shell  electrons.  The  beta  decay  of 
the  fission  fragments  and  the  electrons  released  from  those  decays  are  hypothesized  to  be 
the  main  source  of  the  high-energy  electrons.  This  research  focused  on  the  spatial 
distribution  of  the  fission  fragments  as  the  source  of  high-energy  electrons. 


J.A.Van  Allen,  L.A.  Frank  and  B.J.  O’Brien,  “Satellite  Observations  of  the  Artificial  Radiation  Belt 
of  July  1962,"  Journal  of  Geophysical  Research  68  (1963):  620. 

S  Glasstone  and  P.J.  Dolan.  The  Effects  of  Nuclear  Weapons  (Washington,  DC:  Department  of  the 
Army  DA  Pamphlet  No  50-3,  1977),  27. 
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III.  THE  STARFISH  PRIME  TEST 


A.  PHYSICAL  SETTING 

The  STARFISH  test  oeeurred  on  9  July  1962  at  0900  hours  ZULU  (Greenwieh 
Mean  Time),  whieh  was  2300  hours  in  Hawaii.  The  test  oeeurred  at  16.7°  N  latitude  and 
190.5°  E  longitude  (roughly  above  Johnson  Atoll)  at  an  altitude  of  400  km.38  This  test 
was  eondueted  prior  to  the  implementation  of  the  Limited  Test  Ban  Treaty  (LTBT)  so  at 
that  time,  atmospherie  testing  was  not  prohibited  by  treaty  or  international  law.  The 
announeed  yield  of  STARFISH  was  1.4  MegaTons  (MT).  Most  other  design  details  of 
the  deviee  remain  elassified. 

At  400  km  above  the  Earth’s  surfaee,  the  nighttime  atmosphere  above  Johnston 
Atoll  eonsists  predominantly  of  atomie  oxygen.  NASA  data  for  the  neutral  elements  in 
the  standard  atmosphere  at  nighttime  above  Johnson  Atoll  is  shown  in  Figure  4.  Ionized 
monatomie  oxygen  is  initially  two  orders  of  magnitude  less  than  neutral  oxygen; 
however,  the  prompt  gamma  and  x-ray  energy  released  in  the  detonation  of  a  HANE 
would  ionize  the  eonstituents  of  the  atmosphere  in  the  vieinity  of  the  explosion  (out  to 
10s  of  kilometers).  While  many  gases  eontribute  to  the  atmospherie  eomposition  at  400 
km,  monatomie  oxygen  is  more  abundant  than  other  eonstituents  by  an  order  of 
magnitude.  Therefore,  monatomie  oxygen  is  eonsidered  to  be  representative  of  the 
baekground  gas  for  the  work  in  this  dissertation. 


J.A.  Van  Allen,  L.A.  Frank  and  B.J.  O’Brien,  “Satellite  Observations  of  the  Artificial  Radiation  Belt 
of  July  1962.”  Journal  of  Geophysical  Research  68  (1963):  619. 
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Figure  4.  Number  density  of  gaseous  components  of  nighttime  atmosphere  above 
Johnston  Atoll  as  a  function  of  altitude  (after  Larson,  2011) 

B.  MAGNETIC  BUBBLE  OBSERVATIONS 

The  Earth’s  magnetic  field  is  similar  to  an  offset,  tilted  dipole  field.  It  has  a  dipole 
shape  with  the  southern  pole  of  the  field  occurring  in  the  Earth’s  northern  hemisphere  and 
the  northern  pole  in  the  southern  hemisphere.  Both  magnetic  poles  are  offset  from  the 
respective  geographic  poles  of  the  rotating  Earth.  The  nominal  value  of  the  magnetic  field 
is  0.3  Gauss  at  the  surface  and  decreases  with  distance  further  out  from  the  center  of  the 
Earth’s  surface. 

The  STAREISH  test  employed  27  sensor-equipped  rockets  travelling  at  different 
routes  to  collect  pictures  and  data.  Some  of  these  rockets  collected  data  on  the  magnetic 
field  strength.39  The  data  from  both  particle  counters  and  magnetometers  showed  that  the 
explosion  formed  a  diamagnetic  cavity  around  the  burst  point.  This  asymmetrical  cavity’s 


^9  “A  ‘Quick  Look’  At  the  Technical  Results  of  STARFISH  PRIME,”  6. 
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expansion  and  collapse  occurred  in  about  16  seconds^**  This  magnetic  bubble  is  a 
volume  of  space  where  the  measured  magnetic  field  is  less  than  the  ambient  magnetic 
field.  As  the  Larmor  radius  of  particle  motion  is  a  function  of  magnetic  field  strength,  the 
weaker  magnetic  field  in  the  bubble  causes  a  larger  Larmor,  radius  and  hence  ionized 
particles  travel  in  a  straighter,  less  curved  path  while  traveling  through  the  magnetic 
bubble.  Figure  5  illustrates  the  rocket  trajectories  and  the  time  period  during  which  data 
was  collected  is  shaded  with  red. 


Figure  5.  Time  dependent  size  of  diamagnetic  cavity  produced  by  STARFISFl 

detonation  with  flight  path  of  instrumented  rockets  superimposed  showing 
times  instruments  were  collecting  data  (after  Dyal,  2006) 


C.  OBSERVED  BELT  PUMPING 

The  STARFISH  test  was  centered  on  L=l.  12.^1  The  resulting  magnetic  bubble 
(described  in  the  previous  section)  stretched  the  Earth’s  magnetic  field  lines  and  allowed 


Palmer  Dyal,  “Particle  and  Field  Measurements  of  the  STARFISH  Diagmagnetic  Cavity,”  Journal 
of  Geophysical  Research  111  (2006):  1. 

J.A.  Van  Allen,  L.A.  Frank  and  B.J.  O’Brien,  “Satellite  Observations  of  the  Artificial  Radiation  Belt 
of  July  1962,”  Journal  of  Geophysical  Research  68  (1963):  619. 
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the  ionized  weapon  debris  to  travel  in  straighter  lines  under  kinetic  motion  to  the  edge  of 
the  bubble.  Upon  leaving  the  magnetic  bubble,  the  debris  would  have  re-entered  the 
Earth’s  unperturbed  magnetic  field  and  continued  to  travel  as  specified  under  the  Lorentz 
force  equation  of  motion.  The  magnetic  bubble  would  place  the  ions  a  Larmor  radius 
length  beyond  the  edge  of  the  bubble  where  they  would  resume  their  helical  path  along 
the  magnetic  field  lines  as  previously  described. 

The  United  States  had  several  satellites  with  limited  data  collection  capabilities  in 
orbit  at  the  time  of,  or  immediately  following,  the  STARFISH  test.  These  included  the 
Injun,  Telstar,  Ariel  and  Traac  satellites.  They  all  had  different  instrumentation  and 
orbital  paths,  so  there  is  not  a  common  set  of  data.  Complicating  the  data  comparison,  the 
Telestar  satellite  was  launched  one  day  after  STARFISH  on  10  July  1962.^2  Therefore, 
there  is  no  pre-blast  data  from  Telestar  to  compare  to  the  measurements  taken  post-blast. 

Nevertheless,  all  satellites  reported  an  increase  in  the  flux  of  energetic  electrons 
and  this  flux  deposited  energy  into  the  satellites.  The  energy  from  these  electrons  affected 
the  electronics  and  resulted  in  failures.  Ariel  failed  within  three  days  and  Traac  failed 
within  36  days.^^  An  energetic  electron  impacting  on  a  silicon-based  semiconductor  can 
release  large  numbers  of  secondary  electrons  and  the  corresponding  deposition  of  energy 
can  result  in  damage  to  the  electronic  structure  and/or  memory  failures.  Advances  in 
electronic  technology  make  modem  satellites  even  more  vulnerable  to  upset  or 
destmction  than  technology  from  1962.^4 

The  trapped  particle  data  obtained  underwent  significant  scientific  peer  review 
and  publication.  The  highest  electron  flux  measured  was  1x10^  [electrons/(cm^  sec)]  at 
F=1.20.45  This  peak  flux  fell  to  around  1x10^  [electrons/(cm^  sec)]  at  F=1.8  and  1x10^ 

42  J.A.  Van  Allen,  L.A.  Frank  and  B.J.  O’Brien,  “Satellite  Observations  of  the  Artificial  Radiation  Belt 
of  July  1962,"  Journal  of  Geophyiscal  Research  68  (1963):  622. 

42  Hess,  W.H.,  “The  Artificial  Radiation  Belt  Made  on  July  9,  1962,”  Journal  of  Geophysical 
Research  68  (1963):  667. 

44  Herman  Hoerlin.  United  States  High-Altitude  Test  Experiences,  (Los  Alamos:  LA-6405,  1 976)  3 1 . 

42  Mozer  et  al.,  “Preliminary  Analysis  of  the  Fluxes  and  Spectrums  of  Trapped  Particles  after  the 
Nuclear  Test  of  July  9,  1962”  Journal  of  Geophysical  Research  68  (1963):  648. 
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[electrons/(cm^  sec)]  at  L=2.2 .46  While  this  last  value  is  only  1%  of  the  maximum 
measured  value,  it  is  an  entire  L-shell  beyond  the  location  of  the  detonation;  at  the 
magnetic  equator  that  corresponds  to  a  distance  of  6400  km.  Data  from  the  satellite  Ariel 
indicated  that  there  were  some  high-energy  electrons  out  to  L-shells  5  or  6  but  the  fluxes 
at  those  distances  are  uncertain.47  Explaining  and  modeling  how  the  energetic  electrons 
populated  these  outer  belts  remains  a  scientific  challenge. 


46  W.L. Brown,  W.N.  Hess  and  J.A.  Van  Allen,  “Introduction,”  Journal  of  Geophysical  Research  68 
(1963):  605. 

47  Ibid.  606. 
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IV.  COLLISIONLESS  SHOCK 


A.  FORMATION  OF  COLLISIONLESS  SHOCK 

The  first  hypothesis  eonsidered  was  eollisionless  shoek  aeeeleration  of  the  debris 
ions  as  a  method  to  move  them  to  higher  L-shells  in  the  Earth’s  magnetie  field.  With 
multiple  speeies  in  the  debris,  the  lighter  speeies  eould  drive  the  formation  of  the 
eollisionless  shoek.  As  the  heavier  speeies  passed  through  the  shoek  front,  it  would 
aeeelerate  and  gain  energy  to  move  to  a  higher  L-shell.  After  moving  to  that  higher  L- 
shell,  the  debris  would  undergo  beta  deeay  and  injeet  high-energy  eleetrons  into  the 
Earth’s  magnetie  field  and  eontribute  to  the  data  observed  in  STAREISH. 

Eollowing  a  HANE,  the  ionized  debris  streams  outward  from  the  point  of 
detonation  into  the  Earth’s  magnetie  field.  The  mean  free  path  for  ion-ion  eollisions  is 
large  enough  that  eollisions  between  the  ions  do  not  dominate  the  motion  of  the  partieles 
eompared  to  the  eleetromagnetie  forees.'^^  Collisionless  shoek  is  generated  by  the  highly 
non-linear  interaetions  of  the  expanding  debris  ions,  the  baekground  plasma,  and  the 
eleetromagnetie  fields.  Ultimately  an  eleetrie  field  internets  with  the  baekground  plasma 
and  is  the  meehanism  to  transfer  energy  from  the  debris  ions  to  the  baekground  ions. 
Sinee  the  fields  mediate  the  energy  transfer  instead  of  kinetie  eollisions,  the  proeess  is 
ealled  eollisionless  shoek. 

The  initial  fields  and  early-time  response  of  the  baekground  plasma  ean  be 
analyzed  using  the  eleetron  momentum  equation  with  some  simplifying  assumptions. 
Starting  with  the  eleetron  momentum  equation  as  annotated  by  Hewett,  Earson  and 
Breeht: 

due  —  --  B 

m^n^  —  =  en^E-VP^  +  JeX- 

dt  c  (16) 


R.T.  Glassey  and  W.A.  Strauss,  “Absence  of  Shocks  in  an  Initially  Dilute  Collisionless  Plasma,” 
Communications  in  Mathematical  Physics  113,  (1987):  1. 
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where  nie  is  the  mass  of  the  electron,  He  is  the  number  density  of  the  electrons,  Ue  is  the 
velocity  of  the  electron  fluid,  E  is  the  electric  field,  Pe  is  the  pressure,  Je  is  the  electron 
current  density,  B  is  the  magnetic  field  and  c  is  the  speed  of  light.49  With  the 
approximation  that  electrons  have  zero  mass,  the  equation  reduces  to 

— "  B 

0  =  en^E-yp^  +  JeX- 

c  .  (17) 

In  this  non-relativistic  scenario,  the  zero  inertia  assumption  for  electrons  is  an  appropriate 
assumption  based  on  the  ratio  of  masses  for  electrons  and  nucleons  (electrons  account  for 
0.05%  of  the  simulation  mass).  Re-arranging  the  equation  to  solve  for  E  it  becomes 

-  VP,  JeXB 
en,  en^c 

Then  the  Darwin  approximation  for  Ampere’s  law  is  employed.  The  displacement 
current,  the  time  derivative  of  E  is  ignored.^*'  Writing  the  total  current  density  as  the  sum 
of  the  electron  and  ion  current  densities,  the  Darwin  approximation  to  Ampere’s  law 
becomes 


yxB  =  —J  =  —(je  +  Ji) 
c  c  ^  _ 

Solving  for  Jg  and  substituting  the  results  into  equation  18  gives 


C  /  ^ 

^Pe 

—  Vxfi  -7, 

xB 

eUe  enf 


(19) 


(20) 


D.W  Hewett,  D.J.  Larson  and  S.  Brecht.  “The  Physics  of  Ion  Decoupling  in  Magnetized  Plasma 
Explosions,”  report  no.  LLNL-CONF-  470416  (Livermore,  CA,  Lawrence  Livermore  National  Laboratory, 
2011),  9. 

D.W  Hewett,  D.J.  Larson  and  S.  Brecht.  “The  Physics  of  Ion  Decoupling  in  Magnetized  Plasma 
Explosions,”  report  no.  LLNL-CONF-  470416  (Livermore,  CA,  Lawrence  Livermore  National  Laboratory, 
2011),  9. 
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We  assume  that  early  in  time,  the  gradient  of  the  eleetron  pressure  is  zero  and  the 
magnetie  field  is  irrotational  (has  no  curl)  such  that 


=  0  =  V  X  5 


(21) 


The  quasi-neutral  approximation  applies  when  the  Debye  length  is  small  compared  to  the 
phenomenon  under  investigation.  This  approximation  means  the  number  of  electrons  is 
equal  to  the  number  density  of  the  ions  multiplied  by  their  charge  state^i 


en. 


n  Z  e 

species  species 


species 


(22) 


This  approximation  allows  the  denominator  of  equation  20  to  be  replaced  with  a  sum  of 
the  ion  charge  density,  which  consists  of  debris  ions  and  background  ions.  Also,  early  in 
time  assume  that  the  background  ions  are  stationary 


u 


Br 


(23) 


Now  look  at  the  theta  component  of  the  electric  field  with  the  ambient  magnetic  field  in 
the  z  direction  noting  that  the  subscript  D  refers  to  the  debris  and  the  subscript  B  refers  to 
the  background. 


J  Ir^Z  _  _  j^D^D^Dr)  ^z 

en^c  {Z^n^  +  Z^n,j)c  {Z^n^  +  Z^n,j)c 


(24) 


The  corresponding  equation  for  the  radial  component  of  the  electric  field  is 

— - — - — - 

en^c  (Zj^ns  +  Z^ni^)c  iZ^n,^  +  Z^njj)c  _  ^25) 

Collisionless  shock  is  generated  if  the  angular  component  of  the  electric  field  is 
great  enough  to  accelerate  the  background  ions  ahead  of  the  debris  ions.  The  simple  force 
equation  for  the  acceleration  of  the  background  ions  is 


Francis  Chen,  Introduction  to  Plasma  Physics  and  Controlled  Fusion  Vol  1  (New  York:  Springer, 
2006),  10. 
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(26) 


If  this  ratio  is  greater  than  one,  then  the  eleetrie  field  ean  transfer  enough  energy  to  the 
baekground  to  rapidly  aeeelerate  it  ahead  of  the  debris  and  generate  a  shoek  front  without 
partiele  collisions.  The  effect  is  strictly  electromagnetic. 

In  the  author’s  simulations  of  collisionless  shock,  the  external  magnetic  field  was 
set  in  the  z  direction  with  a  magnitude  of  0.3  gauss.  The  selection  of  the  magnetic  field  in 
the  z  direction  establishes  the  optimal  geometry  with  the  debris  streaming  in  the  radial 
direction  perpendicular  to  the  ambient  magnetic  field.  This  ambient  field  was  affected  by 
the  collisionless  shock  driven  by  the  expanding  plasma.  If  expanding  plasma  generates  a 
collisionless  shock,  it  produces  large  variations  in  the  field  as  seen  in  Figure  6,  which 
shows  the  disturbance  of  the  ambient  magnetic  field  at  a  simulated  time  of  0.08625 
seconds  (or  1725  time  steps)  following  a  HANE. 
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Figure  6.  Collisionless  shoek  observed  in  z  eomponent  of  magnetie  field  during 
author’s  simulation  of  expanding  plasma  pushing  against  baekground 
plasma  after  1725  time  steps  while  the  debris  has  expanded  out  to  10^  cm. 


In  Figure  6,  the  peak  magnetic  field  value  is  0.69  gauss,  which  is  more  than 
double  the  undisturbed  ambient  value.  Also,  the  figure  shows  the  sharp  rise  and  fall  of 
this  peak.  While  this  is  a  snapshot  of  the  magnetic  field,  a  charged  particle  in  space  would 
see  the  time  varying  magnetic  field  in  the  z  direction  and  feel  the  corresponding  electric 
forces  in  the  radial  and  azimuthal  directions. 

B.  COUPLING  OF  DEBRIS  TO  BACKGROUND  PLASMA 

Early  work  on  expanding  plasmas  typically  used  MHD  models  because 
computational  power  was  lacking  for  a  kinetic  ion  model.  MHD  models  rely  heavily  on 
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concepts  from  fluid  dynamics  that  assume  zero  mean-free-path  for  the  ions.  In  a  MHD 
model,  all  the  expanding  plasma  is  assumed  to  eollide  with  the  background  plasma  and 
would  be  contained  by  (or  coupled  to)  the  background  plasma.  The  snowplow  analogy  is 
used  where  the  heavy  snowplow  pushes  against  the  stationary  snow  and  eauses  the  snow 
to  pile  up  in  front  of  the  snowplow;  the  same  effeet  is  seen  in  the  expanding  plasma  with 
the  MHD  model.  The  opposite  analogy  to  the  snowplow  is  the  ship’s  bow.  A  ship  moving 
forward  in  the  water  has  a  narrow  bow  that  pushes  the  water  around  the  ship’s  hull.  This 
motion  creates  a  wake  (or  a  wave)  that  carries  away  some  of  the  energy  but  the  ship 
keeps  moving  forward  and  the  water  does  not  build  up  in  front  of  the  ship.  Those  two 
analogies  ean  help  understand  coupled  (snowplow)  and  uncoupled  (ship’s  bow)  plasmas. 


In  2011,  Hewett  et  al.  demonstrated,  using  a  kinetic  ion  model  for  simulating 
expanding  ionized  debris  into  a  plasma  background,  that  it  was  possible  for  the  debris 
plasma  to  decouple  from  the  background  plasma  and  slip  through  the  background 
plasma.52  They  established  a  eoupling  criterion  given  by 


(^dc  = 


2  y 


(30) 


where  adc  is  the  decoupling  constant,  Z  is  the  charge  state  of  the  debris  (subscript  D)  and 
background  (subscript  E)  and  n  is  the  number  density  of  the  debris  and  background.  If  ttdc 
is  less  than  one,  then  the  debris  will  deeouple  from  the  baekground.  This  slip-through  will 
never  be  captured  by  MHD  models  and  is  potentially  important  beeause  it  is  a  meehanism 
that  allows  debris  ions  to  expand  out  to  approximately  one  gyro-radius  from  the  edge  of 
the  magnetic  bubble. 

The  eoupling  is  best  observed  in  phase  space  plots  of  radial  veloeity  versus  radius. 
Absent  any  baekground  plasma,  the  eharged  partiele  would  spiral  around  a  magnetic  field 
line.  In  phase  space,  this  motion  would  be  represented  by  a  closed  path.  In  the 
simulations  performed  for  this  research  (and  building  off  the  work  done  by  Hewett  et  al), 

^^D.W  Hewett,  D.J.  Larson  and  S.  Brecht.  “The  Physics  of  Ion  Decoupling  in  Magnetized  Plasma 
Explosions,”  report  no.  LLNL-CONF-  470416  (Livermore,  CA,  Lawrence  Livermore  National  Laboratory, 
2011),  9. 
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the  debris  starts  at  the  eenter  of  the  simulation  and  pushes  out  against  the  baekground 
plasma.  In  the  simulations,  the  debris  begins  with  a  radial  veloeity  of  2x10^  cm/s  and  no 
velocity  in  the  z  or  theta  components  such  that 

=  2xl0*[cm/5]  =  _  (31) 

At  the  start  of  the  simulations,  the  debris  is  centered  along  the  axis  (in  cylindrical 
geometry)  at  a  user  specified  radius.  Beyond  the  debris,  the  background  plasma  begins 
and  fills  the  remaining  spatial  domain  of  the  simulation  at  the  user  specified  number 
density.  In  these  simulations,  the  background  plasma  began  with  no  velocity  beyond  its 
thermal  motion.  At  the  end  of  the  simulations,  the  debris  has  slipped  through  of  the 
background  plasma  while  transferring  little  energy  to  the  background  plasma.  Figure  7 
shows  a  typical  simulation  result  in  phase  space  in  which  a  small  increase  in  the  radial 
velocity  of  the  background  plasma  (blue  and  comprised  of  ions)  from  the  original 
value  of  zero  is  observed.  The  radial  velocity  of  the  debris  is  decreasing  because  the 
debris  ions  (represented  by  as  red  and  comprised  of  ions)  are  turning  in  the  ambient 
magnetic  field  and  the  uranium’s  velocity  is  transferred  from  the  radial  component  to  the 
angular  component. 
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Figure  7. 
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Phase  spaee  plot  of  radial  velocity  versus  radius  for  author’s  simulation  of 
uncoupled  plasma  (red)  expanding  into  a  background  plasma  (blue). 


In  the  case  where  the  debris  is  coupled  (by  changing  the  background  density  as 
seen  in  equation  30)  at  the  end  of  the  simulation,  the  debris  (red)  is  seen  to  remain  behind 
the  background  plasma  (blue),  as  shown  in  Figure  8  below.  Additionally,  the  radial 
velocity  of  the  debris  plasma  is  close  to  zero  indicating  that  the  debris  plasma  has  stopped 
expanding  into  the  background  plasma.  Comparing  Figures  7  and  8,  the  effects  of  the 
coupling  can  be  seen  in  the  background  plasma’s  radial  velocity.  In  Figure  8,  the  debris 
plasma  (red)  has  transferred  more  of  its  energy  to  the  radial  velocity  of  the  background 
plasma  (blue)  and  the  background  plasma  is  now  expanding  radial  outward  in  the 
simulation.  The  simulation  involves  no  collision  effects  so  the  transfer  of  energy  from  the 
debris  plasma  to  the  background  plasma  is  purely  electromagnetic. 
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Figure  8.  Phase  spaee  plot  of  radial  veloeity  versus  radius  for  author’s  simulation  of  a 
coupled  debris  plasma  (red)  expanding  into  background  plasma  (blue). 

The  coupled  simulation  generates  results  much  like  those  expected  from  a  MHD 
simulation.53  All  the  energy  from  the  debris  is  expended  in  pushing  the  background 
plasma  and  driving  the  magnetic  field.  The  background  plasma’s  radial  velocity  is  an 
order  of  magnitude  higher  in  the  coupled  simulation  than  the  uncoupled  simulation.  The 
work  by  Hewett  et  al.  provides  a  quick  check  on  whether  the  correct  simulation  tool  is 
MHD  or  PIC.  For  this  work,  PIC  is  the  correct  tool  since  kinetic  effects  are  essential  to 
the  research.  Shocks  with  multiple  species  were  examined  under  the  hypothesis  that  one 
species  would  drive  the  collisionless  shock  and  the  stronger  magnetic  field  generated 
would  couple  species  with  larger  gyroradii. 


Hewett,  D.  Private  conversation  with  author,  23  August  2011. 
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V.  CHARGE  EXCHANGE 


A.  OVERVIEW 

The  process  of  charge  exchange  was  the  second  hypothesis  considered  in  this 
dissertation  research.  Charge  exchange  could  contribute  to  the  motion  of  the  weapon 
debris  across  magnetic  field  lines  and  contribute  to  the  presence  of  the  high-energy 
electrons  at  higher  L-shells  as  observed  following  STARFISH. 

1.  General  Overview 

Charge  exchange  is  the  process  whereby  two  atoms  or  ions  exchange  one  or  more 
electrons  resulting  in  both  changing  their  charge  state.  Because  the  probability  of  charge 
exchange  between  two  ions  is  several  orders  of  magnitude  less  than  the  probability  of 
exchange  between  an  ion  and  a  neutral  atom,  this  work  addresses  only  ion  to  neutral  atom 
charge  exchange.  Also,  the  probability  of  more  than  one  electron  being  exchanged  is 
lower  than  that  of  a  single  electron  exchange^^  and  is  therefore  not  addressed. 
Additionally,  radiative  charge  exchange  is  not  addressed  (in  keeping  with  the  Darwin 
approximation  discussed  in  Chapter  IV.A)  as  the  cross  section  for  radiative  charge 
exchange  is  very  small  compared  to  that  of  other  methods  of  charge  exchange. 

Collisional  charge  exchange  is  usually  analyzed  by  assuming  a  moving  ion  and  a 
stationary  neutral  atom.  These  are  called  the  projectile  and  target,  respectively,  and  that 
nomenclature  is  used  here. 

There  are  two  components  of  charge  exchange  -  electron  capture  and  electron 
loss.  Electron  capture  takes  place  when  the  projectile  captures  an  electron  in  a  collision, 
decreasing  its  charge  state  by  one.  This  has  the  form 

+  (32) 


R.E.  Johnson,  Introduction  to  Atomic  and  Molecular  Collisions  (NY:  Plenum  Press,  1982),  326. 

B.H.  Bransden  and  M.R.C.  McDowell,  Charge  Exchange  and  the  Theory  of  Ion-Atom  Collisions, 
(Oxford:  Clarendon  Press,  1992),  128. 
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where  i  is  the  initial  charge  state  of  the  projectile,  P  refers  to  the  projectile  ion  and  T 
refers  to  the  target  atom.  Electron  loss  occurs  when  the  projectile  loses  an  electron  in  a 
collision,  increasing  its  charge  state.  This  has  the  form 

+  (33) 

where  i  is  the  charge  state  and  /  >  0. 

Physically,  charge  exchange  is  a  quantum  mechanical  process  that  occurs  when 
the  electron  wave  function  overlaps  the  projectile  and  the  target  nuclei.  Experimentally, 
the  cross  section  for  charge  exchange  is  measured  by  shooting  a  beam  of  ions  into  a  gas 
of  known  pressure.  The  ion  beam  is  charge  analyzed  on  the  opposite  side  of  the  gas  and  a 
cross  section  is  calculated.^®  The  probability  of  a  projectile  undergoing  charge  exchange 
is 


p  =  o{E)\v\n^t 

where  p  is  the  probability,  o  is  the  cross  section,  I  v  I  is  the  speed  of  the  projectile  and  n 
is  the  number  density  of  the  target.  The  cross  section  has  the  units  of  area  (cm^  in  these 
calculations)  and  is  energy  dependent.®^  The  easiest  way  to  determine  the  cross  section  is 
experimentally.®^  In  this  work,  the  cross  section  is  calculated  from  the  semi-empirical 
formula  described  below. 


The  cross  section  is  a  method  of  quantifying  the  probability  of  a  specific 
interaction.  In  this  work,  the  interactions  of  interest  are  electron  capture  and  electron  loss, 
and  the  cross  sections  for  those  interactions  are  employed.  The  cross  section  represents  a 
probability  of  interaction  with  a  large  cross  section  representing  a  higher  likelihood  of 
that  specific  interaction. 


®®  H.Knudsen,  H.K.  Haugen  and  P.  Hvelplund,  “Single-electron-capture  cross  section  for  medium-  and 
high-velocity,  highly  charged  ions  colliding  with  atoms”  Physical  Review  A  23  (February  1981):  597. 

R.E.  Johnson,  Introduction  to  Atomic  and  Molecular  Collisions  (HY ■.  Plenum  Press,  1982),  21 

Ibid.  49 
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The  number  density  of  the  target  is  a  funetion  of  radial  distanee.  As  the  distanee 
from  the  eenter  of  the  simulation  inereases,  the  number  density  deereases.  Eventually,  the 
number  density  of  the  target  drops  below  a  threshold  and  the  probability  of  eharge 
exehange  drops  to  a  negligible  amount.  At  this  radial  distanee,  the  eharge  state  is  stable 
and  “frozen  in.”  Charge  state  freeze-out  point  depends  on  the  speed  and  profile  of  the 
number  density  of  the  target  atoms. 

2.  Electron  Capture 

Eleetron  eapture  oeeurs  when  the  projeetile  ion  eaptures  an  eleetron  from  the 
target,  thereby  deereasing  the  eharge  state  by  one.  Eor  an  ion  at  charge  state  +1,  electron 
capture  results  in  the  formation  of  a  neutral  atom.  A  neutral  atom  will  continue  on  its 
trajectory  with  no  influence  from  the  electric  and  magnetic  fields.  If  the  initial  ion  is  at  a 
higher  charge  state  than  +1,  electron  capture  reduces  the  charge  state  by  one  integer.  The 
smaller  charge  state  increases  the  Earmor  radius  and  decreases  the  magnitude  of  the 
effect  of  the  magnetic  and  electric  fields. 

An  ion  capturing  multiple  electrons  in  a  single  event  is  possible  but  less  probable. 
The  probability  of  an  ion  capturing  multiple  electrons  is  generally  an  order  of  magnitude 
lower  than  the  probability  of  a  single  electron  capture. Multiple  electron  captures  in  a 
single  time  step  are  not  simulated  in  this  work. 

3.  Electron  Loss 

Electron  loss  occurs  when  the  projectile  loses  an  electron  to  the  neutral  target  or 
the  continuum,  thereby  increasing  the  charge  state  by  one.  The  important  change  is  to  the 
projectile  charge  state.  A  neutral  atom  that  undergoes  electron  loss  becomes  an  ion  of 
charge  state  +1.  Now  the  motion  of  that  charged  ion  is  affected  by  the  electric  and 
magnetic  fields. 


Hans-Dieter  Betz,  “Charge  States  and  Charge-Changing  Cross  Sections  of  Fast  Heavy  Ions 
Penetrating  Through  Gaseous  and  Solid  Media,”  Review  of  Modern  Physics  44  (1972):  467. 
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A  neutral  atom  losing  multiple  eleetrons  in  a  single  event  is  possible  but  less 
probable.  The  probability  of  a  neutral  atom  losing  multiple  eleetrons  is  generally  an  order 
of  magnitude  less  than  the  probability  of  a  single  eleetron  loss.^*'  Multiple  electron  losses 
in  a  single  time  step  are  not  simulated  in  this  work. 

B,  SEMI-EMPIRICAL  FORMULATION 

The  cross  section  for  charge  exchange  interactions  depends  on  several  variables 
including  the  energy  and  the  charge  state  of  the  projectile.  Cross  sections  can  be 
calculated  using  complex  and  time  consuming  simulations  involving  quantum  mechanics, 
or  they  can  be  measured  experimentally.  A  third  and  simpler  alternative  method  for 
computer  simulations,  especially  appropriate  for  multiple  combinations  of  projectiles  and 
targets  with  a  range  of  energies  and  charge  states,  has  also  been  developed.  Schlachter  et 
al.  originally  published  a  semi-empirical  scaling  law  for  electron  capture  with 
dependencies  on  energy  (velocity),  charge  state  of  the  projectile  and  atomic  number  of 
the  target.61  Schlachter’s  original  scaling  law  has  the  form 

0.5  ^  _1*  10-8  r  .2.2  1  r 

=  m —  l-exp(-0.037*E  )  *  1-exp 

£  L  J  . 

where  q  is  the  charge  state  of  the  projectile,  E  is  the  scaled  kinetic  energy  that  is  a 
function  of  kinetic  energy  per  nucleon  scaled  by  the  atomic  number  and  the  charge  state, 
and  Ztgt  is  the  atomic  number  of  the  target.  Schlachter  states  that  if  q  <  3,  then  the  scaling 
holds  using 

q'  =  q  +  0A 


W.N.Spjeldvik  and  T.A.  Fritz,  “Theory  for  Charge  States  of  Energetic  Oxygen  Ions  in  the  Earth’s 
Radiation  Belts,”  Journal  of  Geophysical  Research  83  (1978):  1585. 

A.S.  Schlachter,  J.W.  Steams,  W.G.  Graham,  K.H.  Berkner,  R.V.  Pyle  and  J.A.  Tanis,  “Electron 
capture  for  fast  highly  charged  ions  in  gas  targets:  An  empirical  scaling  rule,’’  Physical  Review  A  27 
(1983):  3373. 
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as  an  adjustment  to  the  eharge  state  Sehlaehter  et  al.  arrived  at  this  semi-empirieal 
formulation  by  data-fitting  their  experimentally  measured  eross  seetions  to  variables  they 
adjusted  in  their  experiments. 

Shevelko  et  al.  expanded  upon  Sehlaehter’ s  original  work  on  eleetron  eapture  and 
developed  a  semi-empirieal  sealing  law  for  eleetron  loss  as  well.  Their  scaling  law  has 
the  form 


0.5 


tgt 


'^2  U 

l  +  .Ol^  P 

'  u^  +  3.5 

1  31 

4-1-^ — ln(4M-l-l) 


(37) 


where  u  is  the  scaled  velocity  in  atomic  units  (such  that  speed  of  light  corresponds  to  the 
velocity  value  137),  q  is  the  charge  state  of  the  projectile,  Ho  is  the  principle  quantum 
number  of  the  outer  electron  for  the  projectile,  Ztgt  is  the  atomic  number  of  the  target,  and 
Ip  is  the  ionization  potential  for  the  outer  electron  for  the  projectile  atom.^^ 

Using  such  a  semi-empirieal  formula  for  the  cross  section  allows  for  computer 
simulations  across  a  range  of  energies,  charge  states  and  elements.  The  cross  sections  are 
easily  calculated  from  the  known  velocities,  masses  and  charge  states.  This  allows  for 
simulation  of  a  mult-species  (elements)  plasma  with  charge  exchange.  In  a  PIC  code  such 
as  ZMR,  all  the  required  variables  are  known  and  can  be  included. 


C.  MATCH  OF  SEMI-EMPERICAL  FORMULA  TO  EXPERIMENTAL 
DATA 

Sehlaehter  and  Shevelko  do  not  describe  the  exact  method  of  determining  the 
semi-empirieal  fit  in  the  development  of  their  semi-empirieal  formulas.  They  state  that 
their  formulas  were  matched  to  available  experimental  data.  The  lack  of  recent  and 
relevant  data  made  comparing  this  data  to  the  elements  of  interest  difficult.  Recent  data 
for  the  elements  of  interest  at  the  appropriate  velocities  and  energies  was  difficult  to  find. 

62  Ibid.,  3373. 

63  Ibid.,.  3372. 

6^^  V.P.  Shevelko  et  al,  “Charge-changing  processes  in  collisions  of  heavy  many-electron  ions  with  neutral 
atoms,”  Nuclear  Instruments  and  Methods  in  Physics  Research  B  (201 1):  269. 
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Data  from  the  Air  Force  Weapons  Laboratory  (AFWL)  in  1968  covered  many  of  the 
elements  of  interest  at  energies  just  below  the  values  used  in  this  simulation.  The  AFWL 
work  examined  velocities  between  2x10^  cm/s  and  1.2x10*  cm/s.  The  simulations  in  this 
dissertation  research  began  with  an  initial  velocity  of  2x10*  cm/s,  so  the  experimental 
data  set  used  to  benchmark  the  semi-empirical  formulas  does  not  align  exactly  with  the 
energy  spectrum  used  in  the  simulation.  Some  of  the  data  was  reasonably  close  to  the 
semi-empirical  formula  and  some  was  off  by  as  much  as  an  order  of  magnitude.  Figure  9 
is  a  reasonably  close  match  for  electron  capture  of  Xe^^  in  monatomic  oxygen,  which  is 
one  of  the  elements  of  interest  in  this  study. 


Xe  Charge  State  + 1 :  for  electron  capture 
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Figure  9.  Experimental  and  semi-empirical  cross  section  of  electron  capture  of  Xenon 

+1  as  a  function  of  velocity. 
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Xe  Charge  State  + 1 :  for  electron  loss 
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Figure  10.  Experimental  and  semi-empirieal  cross  section  of  electron  loss  of  Xenon  -Hi 

as  a  function  of  velocity. 


Figures  9  and  10  show  a  reasonably  close  correlation  of  experimental  data  to  the 
semi-empirical  formula  for  an  element  of  interest  in  the  correct  energy/velocity  regime  of 
the  simulation.  However,  for  uranium  the  match  between  the  experimental  data  and  the 
semi-empirical  formula  was  too  low  by  almost  an  order  of  magnitude.  Figures  11  and  12 
show  the  plot  of  experimental  data  (points)  against  the  theoretical  cross-section  (line). 
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Cross-section  [cm^] 


Uranium  Charge  State  + 1  for  electron  capture 
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Figure  11.  Experimental  and  semi-empirical  cross  section  of  electron  capture  of 

Uranium+1  as  a  function  of  velocity. 
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Uranium  Charge  State  + 1  for  electron  loss 
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Figure  12.  Experimental  and  semi-empirical  cross  section  of  electron  loss  of 

Uranium+1  as  a  function  of  velocity. 


While  the  velocity  regimes  do  not  perfectly  align,  Schlachter  and  Shevelko’s 
algorithm  is  computationally  efficient.  No  other  approach  to  calculating  cross  sections  for 
both  electron  loss  and  electrons  capture  was  uncovered  in  my  literature  search.  The  code 
was  modified  in  a  modular  manner  and  insertion  of  a  different  algorithm  for  calculating 
cross  sections  would  be  straightforward.  Finding  a  more  accurate  algorithm  for 
calculating  cross-sections  for  electron  capture  and  electron  loss  in  the  regimes  of  interest 
remains  a  goal. 
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VI.  ZMR  CODE  DESCRIPTION 


A.  OVERVIEW 

ZMR  is  a  two  dimensional,  plasma  simulation,  FORTRAN  code  created  by  Dr. 
Dennis  Hewett  of  LLNL.  It  uses  the  particle-in-cell  (PIC)  technique  in  cylindrical 
geometry  to  simulate  plasma  behavior.  It  has  been  used  for  over  30  years  on  topics  that 
include  theta  pinch  plasma,  z  pinch  plasma,  and  HANE  debris  coupling. 

ZMR  is  a  hybrid  code  that  treats  the  ions  in  the  PIC  methodology  and  the 
electrons  as  an  inertia-less  charged  fluid.  As  the  ratio  of  electron  to  proton  mass  is 
1/1836,  less  than  0.05%  of  the  mass  in  the  non-relativistic  simulation  is  electrons,  so  this 
is  a  valid  assumption.  It  is  assumed  that  the  ions  contain  all  of  the  inertia,  and  the 
electrons  are  a  mass-less  fluid  that  flows  along  with  the  heavy  ion  particles.  Hybrid  codes 
provide  good  resolution  for  the  kinetic  effects  of  the  ions,  but  cannot  resolve  the  kinetic 
effects  of  the  electrons.®®  A  full  PIC  code  that  resolves  the  electron  motion  requires  a 
very  small  spatial  domain  to  resolve  the  electron  Debye  length  and  a  very  small  time  step 
to  resolve  the  electron  plasma  frequency.®^  For  plasma  at  STARFISH  altitudes  with  an 
energy  of  1  electron-volt  (eV),  this  resolution  requires  a  cell  length  of  0.07  cm  and  a  time 
step  on  the  order  of  10~"  seconds.  This  level  of  resolution  is  computationally  expensive 
and  six  orders  of  magnitude  greater  than  what  is  required  for  a  hybrid  code.  Since  this 
work  focuses  on  the  spatial  distribution  of  the  debris  ions,  using  a  hybrid  code  provides 
the  appropriate  balance  of  fidelity  in  spatial  domain  and  speed  of  computing. 

The  code  also  employs  the  quasi-neutral  approximation  of  plasma  physics,®^ 
which  is  applicable  when  the  Debye  length  is  much  smaller  than  the  length  scales  under 
investigation.  The  Debye  length  is  a  measure  of  the  length  over  which  the  plasma  will 

®®  D.W.  Hewett.  “Documentation  of  the  EM  2.5  Quasi-neutral  Hybrid  Code  ZEMER,”  p3 

®®  Petr  Hellinger,  “Structure  and  Stationarity  of  Quasi-perpendicular  Shocks:  Numerical  simulations.” 
Planet.  Space  Sci,  51  (2003):  650. 

®7  Ibid.  651. 

®^  D.W.  Hewett.  “Documentation  of  the  EM  2.5  Quasi-neutral  Hybrid  Code  ZEMER,”  p2 
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shield  the  field  introdueed  by  a  test  eharge.  If  this  eondition  holds,  the  plasma  is 
eonsidered  “neutral  enough”  and  mathematieally  the  electron  density  is  approximately 
equal  to  the  ion  density. is  written  as 

.  (38) 

This  approximation  becomes  important  when  calculating  electrical  currents  as  well  as 
solving  the  field  equations. 

The  geometry  assumes  symmetry  in  the  theta  (angular)  direction  such  that  the 
partial  derivative  of  any  variable  with  respect  to  6  is  zero. 

(39) 

and  the  particle  displacement  in  the  theta  direction  is  set  to  zero  in  every  time  step.^*' 
Thus,  the  direction  of  interest  in  ZMR  is  radial.  However,  the  theta  components  of  all  the 
vectors  are  calculated  and  maintained  throughout  the  simulation.  For  this  reason,  ZMR  is 
referred  to  a  2.5  dimension  code. 

Different  problems  are  modeled  in  ZMR  by  different  inputs  at  the  start  of  the 
simulation.  The  input  module  allows  the  user  to  set  the  spatial  parameters  of  the 
simulation  as  well  as  the  size  of  the  time  step.  Within  the  input  module,  the  user  can 
specify  the  number  of  species  modeled  in  the  simulation,  as  well  as  the  starting 
geometries,  densities  and  velocities  of  each  species.  The  beginning  geometry,  density  and 
atomic  mass  of  a  species  determine  the  initial  mass  of  that  individual  species. 

Since  the  total  mass  of  a  species  is  specified  and  the  total  number  of  simulation 
particles  for  that  species  is  also  specified,  ZMR  computes  the  number  of  real  particles  per 
simulation  particle  (RPPSP).  This  quantity  is  initialized  for  the  species  but  is  tracked  for 
each  individual  simulation  particle.  As  a  simulation  particle  quantity,  it  enables  the 


Francis  Chen,  Introduction  to  Plasma  Physics  and  Controlled  Fusion  Vol  1  (New  York:  Springer, 
2006),  10. 

D.W.  Hewett,  Documentation  of  the  EM  2.5  Quasi-neutral  Hybrid  Code  ZEMER,  3. 
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addition  of  the  charge  exchange  routine  by  allowing  a  simulation  particle  to  “split”  into 
two  simulation  particles.  The  simulation  particle  that  splits  from  the  original  simulation 
particle  will  have  the  same  location  and  velocity  but  different  charge  states  reflecting  that 
some  of  the  real  particles  modeled  by  the  simulation  particle  have  undergone  charge 
exchange.  The  total  mass  of  the  simulation  is  unchanged  as  the  original  simulation 
particle  has  its  RPPSP  value  decreased  by  the  same  amount  as  the  new  particle  created 
with  a  different  charge  state.  This  procedure  conserves  mass,  energy  and  momentum 
within  the  simulation  but  it  does  result  in  an  exponential  increase  in  the  number  of 
simulation  particles.  As  the  number  of  simulation  particles  increases,  the  simulation 
undergoes  more  individual  calculations  and  slows  down.  Including  charge  exchanged 
increased  the  simulation  run  time  (wall  time)  by  a  factor  of  ten. 

Advancing  the  position  is  done  using  the  Boris  Push  method.  The  Boris  Push 
method  is  an  algorithm  to  time  advance  the  Lorentz  force  equation.  Generally,  the  Boris 
Push  method  uses  a  half  time  step  advancement  from  the  electric  field. 


=  -  +S-E  — 

^  prime  ^  old 

m  2 


(40) 


a  full  time  step  advancement  from  the  magnetic  field 
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followed  by  a  rotation  of  that  magnetic  field  advancement 
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and  then  the  last  half  time  step  advancement  from  the  electric  field 
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A  more  detailed  explanation  of  the  Boris  Push  method  is  available  in  Birdsall  and 
Langdon.^i  The  new  positions 


^new  ^old 


+  v„„.  At 


(44) 


of  all  the  particles  are  used  to  calculate  the  charge  densities  and,  along  with  the  new 
velocities,  are  used  to  calculate  the  electric  currents. 


After  the  particles  are  advanced  in  time,  the  electric  and  magnetic  fields  are 
advanced  in  time.  The  cylindrical  geometry  allows  the  theta  component  of  the  vector 
potential  to  be  solved  for  explicitly  from  the  electron  momentum  equation  (in  the  limit  of 
small  electron  inertia).  From  the  explicitly  calculated  value  of  the  theta  component  of  the 
vector  potential,  the  radial  and  z  components  of  the  magnetic  field  are  calculated.  The 
theta  component  of  the  magnetic  field  is  calculated  in  a  manner  similar  to  the  theta 
component  of  the  vector  potential.  The  components  of  the  electric  field  are  calculated 
from  the  electron  momentum  in  equation  18. 

ZMR  employs  the  Darwin  approximation  that  removes  the  propagation  of 
electromagnetic  waves. This  approximation  means  changes  to  the  electric  and 
magnetic  fields  propagate  instantly  across  the  simulation  with  no  time  delay.  This  is  a 
valid  assumption  because  the  study  is  focused  on  the  motion  of  ions  and  not  the 
electromagnetic  waves  produced  by  the  plasma.  The  electric  and  magnetic  fields  are 
computed  and  stored  on  the  simulation  grid.  Only  the  closest  four  grid  points  contribute 
to  the  fields  used  in  the  equations  of  motion  for  the  particle  at  a  particular  point.  In 
Figure  13  below,  the  electric  and  magnetic  fields  at  the  simulation  particle  represented  by 
the  green  dot  are  calculated  from  the  values  of  those  fields  stored  on  the  grid  represented 
by  the  blue  dots.  The  other  values  of  the  fields  on  the  grid  (represented  by  the  red  dots) 
do  not  impact  the  motion  of  simulation  particle.  In  the  computational  time  step,  an 


C.K.  Birdsall  and  A.B.  Langdon,  Plasma  Physics  Via  Computer  Simulation.  (New  York:  Taylor  and 
Francis  Group,  2005),  59. 

M.R.  Gibbons  and  D.W.  Hewett,  “The  Darwin  Direct  Implicit  Particle-in-Cell  Method  for 
Simulation  of  Low  Frequency  Plasma  Phenomena,”  Journal  of  Computational  PhysicsllO  (1995):  232. 
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electromagnetic  wave  would  prorogate  across  four  cells  so  assuming  the  wave  is  faster 
does  not  change  the  physics. 


Figure  13.  View  of  simulation  grid  with  one  simulation  particle  (green)  and  diamonds 
showing  cell  centered  field  values.  Only  the  blue  diamonds  contribute  to  the 
force  computations  for  the  green  simulation  particle. 


The  data  from  ZMR  is  captured  in  output  data  files.  For  every  species,  the 
position  and  velocities  are  recorded.  The  magnetic  field  is  also  captured  to  observe  the 
magnetic  shock  front. 


B,  COMPUTATIONAL  FLOW  IN  ZMR 

ZMR  follows  the  general  flow  of  hybrid  PIC  codes  of  advancing  the  particles 
using  the  Lorentz  force,  accumulating  the  source  terms  on  the  grid,  calculating  the 
electric  and  magnetic  fields  and  then  advancing  the  time  step.  The  following  is  a  brief 
overview  of  all  those  steps. 

At  the  beginning  of  a  time  step,  the  particles  are  moved  according  to  the  Lorentz 

force 


dv 

dt 


V  new  V  old 

At 


=  q(E  +  vxB) 


(45) 


This  vector  equation  can  is  decomposed  into  its  scalar,  directional  components  such  as 
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for  the  X  component  and  similar  expansions  for  the  y  and  z  components.  This  difference 
equation  is  solved  for  the  time-advanced  velocity  in  all  three  components.  Then  the 
positions  of  all  the  particles  are  advanced.  From  this  time  advanced  positions,  the 
currents  and  charge  densities  are  calculated  and  assigned  to  the  internal  grid  structure. 

After  the  time-advanced  positions,  currents  and  charge  densities  are  calculated, 
the  time-advanced  values  for  the  electric  and  magnetic  fields  are  calculated.  ZMR 
capitalizes  on  the  ability  to  solve  for  the  theta  component  of  the  vector  potential  (A0) 
explicitly.  The  angular  symmetry  required  in  ZMR  and  the  use  of  cylindrical  geometry 
allows  for  the  closed  form  solution  to  Ag.  Once  the  theta  component  of  the  vector 
potential  (A)  is  found,  the  radial  and  axial  components  of  A  are  calculated.  From  the 
components  of  A,  the  components  of  the  magnetic  field  (B)  are  easily  calculated  using 

V  X  A  =  B  (47) 

since  all  the  components  of  the  vector  potential  are  known.  A  complete  explanation  of  the 
field-solving  algorithm  is  available  from  Hewett.^3  computer  program  iterates 

position  advancement  and  then  a  field  calculation  every  time  step.  The  number  of  time 
steps  is  a  user  specified  input. 

The  species  of  interest  in  this  work  were  derived  from  selected  debris  and  fission 
fragments  from  the  nuclear  weapon  and  the  ambient  ionosphere.  As  shown  in  Figure  4, 
atomic  oxygen  is  the  dominant  constituent  of  the  unperturbed  atmosphere  at  the  time  and 
elevation  of  STARFISH,  so  singly  ionized  oxygen  was  used  as  the  background  plasma 
and  neutral  oxygen  was  used  as  the  catalyst  for  charge  exchange.  For  fission  fragments, 
the  elements  of  xenon  and  zirconium  were  used  to  correspond  to  the  high  points  of  the 
fission  fragment  distribution  curve.  Lead  and  carbon  were  also  examined  as  elements  at 
the  extreme  heavy  and  light  ends  of  the  mass  spectrum. 

D.W.  Hewett,  “A  Global  Method  of  Solving  the  Electron-Field  Equations  in  a  Zero-Inertia- 
Electron-Hybrid  Plasma  Simulation  Code,”  Journal  of  Computational  Physics  38  (1980):  379-395. 
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C.  INITIAL  CONDITIONS  USED  IN  THE  SIMULATIONS 

In  this  work,  a  common  set  of  user-supplied  physical  parameters  was  used  in 
order  to  allow  a  eonsistent  comparison  of  results.  The  grid  was  divided  into  256  cells  in 
the  radial  direetion  and  one  in  the  axial  direction  (plus  two  ghost  cells  in  the  axial 
direetion  to  calculate  the  fields  primary  cell).  The  spatial  domain  in  the  radial  direction 
was  set  to  600  km.  The  spatial  domain  in  the  axial  direction  was  set  to  10  km.  The  time 
step  was  set  to  5x10“^  seeonds.  These  parameters  were  used  throughout  all  the 
simulations  examining  both  collisionless  shoek  generated  from  multiple  speeies  and 
eharge  exehange. 

D,  SPECIFIC  MODIFICATIONS  FOR  THIS  WORK 

The  ZMR  Code  was  the  starting  point  for  this  work.  The  simple  geometry  and 
speed  of  execution  made  it  a  natural  starting  point  to  examine  the  physics  of  collisionless 
shoek  and  eharge  exehange.  Sinee  eharge  exchange  appears  to  have  a  significant  impact 
on  the  motion  and  final  spatial  distribution  of  the  debris,  it  is  being  incorporated  into 
other  codes  at  LLNL  or  other  researeh  facilities  that  have  higher  fidelity  and  more 
representative  geometry. 

In  this  work,  the  intent  was  to  use  ZMR  to  simulate  a  1-D  radial  slice  of  the 
expanding  plasma.  In  order  to  reduce  the  spatial  dimensionality,  the  simulations  were 
initialized  with  only  one  cell  in  the  z  direction.  That  one  cell  width  was  set  to  10  km, 
large  enough  to  prevent  simulation  particles  from  leaving  the  simulation  through  the  ends 
of  the  cylindrical  mesh.  These  dimensions  are  large  enough  to  replicate  the  distance  of 
relevance  for  STARFISH  and  simulate  the  motion  of  the  debris. 

1.  Modifications  for  Including  Multiple  Species 

ZMR  was  run  using  a  debris  plasma  eonsisting  of  multiple  speeies.  This  ability 
was  already  resident  within  ZMR  but  Hewett’s  previous  work  foeused  on  single  speeies 
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debris  for  his  studies  of  eoupling^^  Using  multiple  debris  speeies  required  only  slight 
modifieations  to  the  input  file. 

In  the  initial  simulations,  the  two  speeies  were  different  eharge  states  of  uranium. 
One  speeies  maintained  a  eharge  state  of  +1  while  the  seeond  speeies  had  a  charge  state 
of  +2,  +4  or  +8,  depending  on  the  simulation.  The  ratio  of  charge  state  to  mass 

determines  the  Larmor  radius  (as  defined  in  equation  14).  Holding  the  mass  constant 
while  varying  the  charge  state  allowed  for  different  plasma  interactions  while  simplifying 
the  number  of  variables.  Keeping  a  constant  mass  and  number  density  in  the  simulation 
also  kept  the  kinetic  energy  constant.  It  was  possible  to  understand  the  effects  of  the 
multiple  species  on  the  formation  of  collisionless  shock  by  varying  the  charge  state  of  a 
single  element/isotope. 

The  debris  plasmas  were  initialized  in  two  different  configurations.  In  the  first 
configuration,  all  the  debris  plasma  was  co-located  in  a  single  mix  with  a  radius  of  10 
km.  In  the  simulation  the  debris  plasmas  expanded  into  the  background  plasma  and  the 
ambient  magnetic  field  and  the  different  Larmor  radii  separated  the  two  charge  states.  In 
the  second  configuration,  the  debris  plasma  was  initially  stacked  into  two  shells  with  the 
higher  charge  state  on  the  outside  with  an  outer  radius  of  10  km.  The  lower  charge  state 
was  on  the  inside  of  the  expanding  cylinder  with  an  initial  radius  of  5km.  This  geometry 
relates  to  the  heavier  element  on  the  inside  and  the  lighter  element  on  the  outside  if  all  the 
charge  states  were  set  to  +1.  Again  the  debris  species  expanded  into  the  background 
plasma  and  the  ambient  magnetic  field. 

All  of  the  multispecies  simulations  without  charge  exchange  were  conducted 
using  a  constant  atmosphere.  The  constant  density  was  used  to  focus  on  the  effects  from 
multiple  species  on  the  collisionless  shock  formation  and  not  on  the  atmospheric  density. 

2,  Modifications  to  Analyze  Charge  Exchange 

ZMR  had  to  be  modified  to  include  the  physics  of  charge  exchange.  The  thrust  of 
the  research  was  to  determine  if  charge  exchange  would  affect  the  final  spatial 

D.W.  Hewett,  private  conversation  with  author,  28  July  2011. 
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distribution  of  the  debris.  Could  neutral  debris  produced  by  charge  exchange  make  it 
through  the  shock  front?  If  so,  how  much  and  at  what  velocity?  The  desire  to  investigate 
these  effects  guided  how  the  physics  of  charge  exchange  was  implemented  in  the  code 
modifications. 

The  two  competing  processes  of  charge  exchange  (electron  capture  and  electron 
loss)  could  be  dealt  with  separately  or  combined  into  a  single  process.  Combining  them 
into  a  single  process  would  involve  comparing  the  rates  of  electron  capture  and  electron 
loss  and  calculating  the  difference  between  the  two  rates.  This  difference  would  be 
manifested  as  a  change  in  the  charge  state.  This  method  does  not  reflect  the  two-step 
charge  exchange  process.  First,  a  charged  ion  undergoes  electron  capture  and  becomes 
neutral.  Then  that  neutral  atom  travels  a  distance  before  undergoing  electron  loss  and 
becoming  charged  again.  The  distance  traveled  as  a  neutral  atom  is  a  kinetic  effect  that 
will  determine  the  final  spatial  distribution  of  the  debris  and  thus  could  not  be  discounted. 

Another  option  for  modeling  charge  exchange  would  involve  adjusting  the  charge 
state  of  the  simulation  particle  to  account  for  the  relative  change  in  charge  state.  The 
charge  state  of  the  species  is  an  input  parameter  and  is  typically  an  integer.  However,  in 
ZMR  the  variable  for  charge  state  is  defined  as  a  real  number,  and  thus  could  have  a 
fractional  value.  Based  on  the  difference  in  the  rate  of  charge  exchange,  the  variables  for 
charge  state  could  be  adjusted  to  account  for  the  difference  in  electron  capture  and 
electron  loss  cross  sections  computed.  Adjusting  the  charge  state  would  change  the 
Larmor  radius  for  that  species.  This  adjustment  would  also  affect  the  final  spatial 
distribution;  however  it  would  not  appropriately  simulate  physical  behavior.  Ions  do  not 
have  charge  states  with  fractional  values  and  this  technique  would  produce  charge  states 
with  non-integer  values.  Therefore,  this  computational  approach  would  have  produced 
non-physical  results  and  was  not  pursued. 

The  modification  used  in  this  work  involved  transporting  the  simulation  particles 
through  neutral  oxygen.  Each  simulation  particle  represents  many  real  particles  and  is 
tracked  by  a  variable  assigning  the  number  of  real  particles  per  simulation  particle 
(RPPSP).  A  certain  number  of  these  real  particles  will  undergo  charge  exchange.  That 
number  was  calculated  directly  using 
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number  =  <7{E)\v\ngR(At)  (48) 

where  a(E)  is  the  ealculated  eross  section  for  either  electron  capture  or  electron  loss,  |v|  is 
the  speed,  At  is  the  time  step,  Uo  is  the  number  density  of  neutral  oxygen,  and  R  is  the 
number  of  real  particles  per  simulation  particle.  The  existing  simulation  particle  is 
decremented  in  how  many  simulation  particles  it  represents  (R)  and  a  new  simulation 
particle  of  different  charge  state  is  created  with  the  same  spatial  properties  as  the  original 
simulation  particle.  The  total  number  of  real  particles  in  the  simulation  is  unchanged.  The 
newly  created  simulation  particle  has  an  integer  value  for  charge  state  and  a  clearly 
defined  Larmor  radius;  this  would  not  be  the  case  if  the  charge  states  were  just  adjusted 
to  a  fraction.  This  new  simulation  particle  will  move  at  least  a  time  step  before 
undergoing  charge  exchange  and  this  movement  will  affect  the  spatial  distribution  of  the 
debris.  This  method  of  tracking  charge  exchange  creates  the  problem  of  geometric 
growth  in  the  number  of  simulation  particles.  Use  of  the  high  performance  computing 
resource  at  LLNL  allowed  for  relatively  quick  runs  using  two  debris  species. 

The  charge  exchange  routines  were  inserted  into  the  original  code  in  the  particle 
movement  phase  of  the  logic  flow.  After  all  the  particles  have  their  new  locations 
determined,  the  charge  exchange  routines  are  applied  to  the  debris  ions.  To  control  the 
geometric  increase  in  the  number  of  simulation  particles,  all  debris  particles  start  at 
charge  state  +1.  No  higher  initial  charge  states  were  considered  in  this  work  examining 
the  effects  of  charge  exchange.  Those  ions  at  charge  state  +1  underwent  electron  capture 
and  those  neutral  debris  atoms  underwent  electron  loss.  Once  all  the  charge  states  were 
adjusted,  then  ZMR  collected  the  source  terms  of  charge  density  and  current  to  be  used  in 
the  subsequent  field  calculations.  The  fields  are  calculated  using  the  source  terms 
appropriately  adjusted  for  charge  exchange,  and  reflect  the  correct  values  for  each  time 
step.  The  value  of  the  time  step  (At)  needed  to  be  kept  small  enough  so  that  not  all  the 
particles  underwent  charge  exchange  in  a  single  time  step.  Having  all  the  atoms  undergo 
charge  exchange  in  a  single  time  step  is  not  physically  realistic  and  would  introduce 
instabilities  in  the  calculations. 
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Another  modification  involves  the  use  of  an  exponentially  declining  atmospheric 
density.  The  particle  density  for  the  atmosphere  is  dependent  on  the  height  above  the 
surface  of  the  Earth.  This  decrease  with  altitude  is  exponential  and  occurs  for  both  the 
neutral  and  ionized  oxygen  atoms.  After  modification,  ZMR  could  be  run  with  either  a 
constant  particle  density  or  an  exponentially  decreasing  atmosphere.  The  number  density 
of  the  neutral  oxygen  species  is  used  only  to  compute  the  amount  of  charge  exchange;  the 
neutral  oxygen  species  has  no  other  impact  on  the  plasma  dynamics.  Therefore,  neutral 
oxygen  atoms  are  not  represented  with  simulation  particles  but  its  density  is  stored  on  the 
grid  of  the  simulation  and  called  at  the  appropriate  value  for  the  location  of  the 
simulation  particle  undergoing  charge  exchange.  The  number  density  of  the  neutral 
oxygen  atoms  could  be  made  an  exponentially  decreasing  function  of  radial  distance  to 
simulate  the  exponential  decrease  in  real  density. 

It  was  assumed  that  all  the  initially  neutral  oxygen  within  50  km  of  the  detonation 
was  highly  ionized  by  the  detonation  and  no  neutral  oxygen  existed  within  this  range 
following  the  detonation.  In  both  the  constant  density  and  exponential  atmosphere,  the 
charge  exchange  routines  were  not  called  while  the  simulation  particle’s  displacement 
from  the  centerline  axis  was  less  than  50  km. 

Unlike  neutral  oxygen  atoms,  ionized  oxygen  in  the  background  directly  affects 
the  kinetics  of  the  debris  plasma.  Therefore,  the  ionized  oxygen  is  represented  by 
simulation  particles  and  ZMR  computes  all  the  kinetics  for  ionized  oxygen  as  the 
background  plasma.  The  exponentially  decreasing  number  density  for  ionized  oxygen 
was  simulated  by  exponentially  decreasing  the  real  particles  per  simulation  particle 
(RPPSP)  variable  assigned  to  the  ionized  oxygen  simulation  particles.  Figure  14  shows 
the  number  density  of  the  real  background  plasma  particles  used  in  ZMR  for  the 
exponentially  decreasing  atmosphere.  This  allowed  the  simulation  of  the  interaction 
between  the  debris  plasma  and  the  background  plasma  to  capture  the  kinetic  effects  while 
including  an  exponentially  decreasing  atmosphere. 
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Figure  14.  Oxygen+1  background  plasma  number  density  as  a  function  of  radius. 

E.  VERIFICATION  OF  ZMR 

ZMR  has  been  used  for  the  simulations  of  plasma  dynamics  for  over  30  years  at 
LLNL.  While  it  uses  a  different  algorithm  to  solve  for  the  electric  and  magnetic  fields,  its 
results  are  comparable  to  those  from  other  hybrid  plasma  codes.  During  the  modification 
of  ZMR  in  this  dissertation  research,  care  was  taken  to  verify  that  results  from  reference 
calculations  were  the  same  before  and  after  each  modification.  Additionally,  cross  section 
calculations  from  ZMR  match  cross  section  calculations  performed  by  hand.  This 
provided  a  high  level  of  confidence  that  the  code  modifications  were  implemented 
correctly. 
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As  a  separate  test,  ZMR  was  also  modified  to  remove  the  plasma  background  and 
replace  it  with  only  a  neutral  hydrogen  background.  The  debris  plasma  was  injected  as  a 
beam  directed  into  the  neutral  hydrogen.  According  to  the  Trapped  Radiation  Handbook, 
the  fraction  of  atoms  that  should  be  neutral  is  given  by 


(49) 

where  Y  is  the  neutral  fraction,  ooi  is  the  ionization  cross  section,  oio  is  the  neutralization 
cross  section,  p  is  the  mass  density,  s  is  the  penetration  depth  and  m  is  the  atomic  mass.^^ 
The  steady  state  case  allows  for  very  deep  penetration  of  the  projectile  ion  into  the 
neutral  background  and  s  is  very  large.  The  steady  state  value  of  the  neutral  fraction 
reduces  to  just  the  leading  term  of  equation  49.  To  run  this  test  case  computationally,  the 
external  magnetic  field  was  turned  off  (Bz=0)  so  the  uranium  ions  would  not  turn.  The 
background  plasma  was  set  to  neutral  atomic  hydrogen  so  there  was  no  coupling  and  no 
energy  transfer  between  the  uranium  ions  and  the  background.  The  exponential  neutral 
atmosphere  profile  was  replaced  with  a  constant  density  profile.  At  a  velocity  of  2x10^ 
cm/s,  the  cross  section  for  electron  capture  is  calculated  as  2.7x10“^^  cm^  and  the  cross 
section  for  electron  loss  is  calculated  as  2.75x10”'^  cm^.  The  asymptotic  limit  for 
equation  49  is  0.48.  Figure  15  show  Equation  49  plotted  as  a  mathematical  function. 
Figure  16  shows  the  fraction  of  neutral  particles  from  the  simulation  approaching  the 
theoretical  limit  of  0.48.  The  simulation  data  in  Figure  16  also  has  the  calculated 
characteristic  time  (e'*)  of  0.003  seconds.  The  similarity  between  Figures  15  and  16  show 
the  similarities  between  the  physics  principle  and  the  computational  simulation  of  that 
principle.  This  similarity  furthers  the  confidence  that  modification  to  ZMR  to  include 
charge  exchange  was  implemented  correctly. 
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Figure  15.  Neutral  fraction  of  uranium  in  theoretical  calculation  as  a  function  of 

penetration  distance. 
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Figure  16.  Neutral  fraction  of  uranium  in  test  simulation  as  a  function  of  real  time. 

This  match  between  the  theoretical  calculations  and  the  simulation  results  for  a 
simple  case  gives  higher  confidence  that  the  simulation  is  correctly  applying  the  same 
physics  in  a  more  complicated  case  with  velocity  dependent  cross  sections. 

F.  OVERVIEW  OF  SIMULATIONS  EMPLOYED 

This  work  employed  simulation  of  plasma  to  reveal  the  underlying  physics  of  the 
debris  streaming  outward  from  a  FIANE.  Because  the  principles  under  investigation  are 
kinetic  in  nature  and  involved  the  movement  of  ions  under  the  influence  of  magnetic 
fields,  including  those  that  are  self-induced,  a  particle-in-cell  (PIC)  code  was  employed. 
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This  work  relied  primarily  on  a  previous  eode  developed  at  LLNL  -  ZMR,^^  with 
modifieations  to  introduee  additional  physios.  The  resulting  modifioations  of  ZMR 
provided  the  vehiole  to  examine  two  different  prooesses  -  (1)  oollisionless  shook  driven 
by  multiple  species;  and  (2)  charge  exchange,  as  well  as  the  effects  of  those  processes  on 
the  final  spatial  disposition  of  the  debris  from  a  HANE.  These  two  processes  were 
considered  as  possible  explanations  for  the  anomalous  debris  distributions  observer  in 
tests  such  as  STARFISH. 

The  first  process  investigated  was  collisionless  shock  driven  by  multiple  species 
of  debris  ions.  To  keep  the  comparisons  similar  across  various  simulations,  the 
simulations  considered  multiple  ionic  states  of  the  single  isotope  uranium-238.  The 
effects  of  multiple  species  were  introduced  by  varying  the  charge  state  of  the  uranium 
atoms.  One  debris  species  of  uranium  was  kept  at  charge  state  +1  while  the  second  debris 
species  of  uranium  had  a  variable  charge  state,  taking  values  of  +2,  +4  and  then  +8.  The 
variable  charges  served  as  a  surrogate  for  consideration  of  other  lighter  elements.  The 
advantage  of  varying  charge  state  in  lieu  of  mass  was  consistency  in  comparison  between 
simulations.  Since  the  mass  and  number  density  were  the  same  between  simulations,  the 
various  simulations  could  be  conducted  with  the  same  total  energy  content.  The  different 
gyroradii  (from  the  different  charge  states)  was  the  variable  under  consideration. 

This  evaluation  of  collisionless  shock  simulated  two  different  configurations  of 
the  debris.  In  the  first  configuration,  the  debris  species  were  considered  to  be  mixed 
evenly  at  the  starting  radius  (i.e.,  the  point  of  detonation).  This  first  configuration  was 
simulated  to  serve  as  a  proof  of  principle  and  to  confirm  the  proper  functioning  of  the 
model.  The  ZMR  program  was  programmed  to  handle  two  debris  species  (with  a  set 
number  density)  in  which  the  user  specifies  the  initial  region  of  the  debris.  In  the  second 
configuration,  the  debris  was  staggered  in  initial  distribution;  the  outer  debris  had  the 
lower  charge  state  (charge  state  +1),  representing  the  heavier  ions.  The  inner  debris  had  a 
higher  charge  state  (charge  state  +2,  +4  or  +8),  representing  the  lighter  ions.  These 
second  configuration  simulations,  with  the  staggered  debris,  tested  the  hypothesis  that  the 

D.W.  Hewett,  Documentation  of  the  EM  2.5  Quasi-neutral  Hybrid  Code  ZEMER,  1. 
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large  gyroradii  ions  could  drive  collisionless  shock  disproportionally  affecting  the  smaller 
gyroradh  (i.e.,  lighter  or  more  highly  charged)  ions  as  they  passed  through  the  shock. 
Table  1  summarizes  the  simulations  carried  out  to  evaluate  the  hypothesis  of  collisionless 
shock. 


Simulation 

Debris 
Species  1 

Debris 
Species  2 

Debris 

Configuration 

Hypothesis 

tested 

CS-1 

U+1 

U+1 

Co-located  - 
evenly  mixed 

Proof  of  coding 

CS-2 

U+1 

U+2 

Co-located  - 
evenly  mixed 

Effects  of  lighter 
species  -  half 
mass 

CS-3 

U+1 

U+4 

Co-located  - 
evenly  mixed 

Effects  of  lighter 
species  -  quarter 
mass 

CS-4 

u+1 

U+8 

Co-located  - 
evenly  mixed 

Effects  of  lighter 
species  -  eighth 
mass 

CS-5 

u+1 

U+2 

Staggered 

Effects  of  lighter 
species  -  half 
mass 

CS-6 

u+1 

U+4 

Staggered 

Effects  of  lighter 
species  -  quarter 
mass 

CS-7 

u+1 

U+8 

Staggered 

Effects  of  lighter 
species  -  eighth 
mass 

Table  1.  List  of  simulations  conducted  to  examine  collision  shock  driven  by  multiple 

species. 


Following  the  simulations  with  multiple  species  focusing  on  the  collisionless 
shock  mechanism,  the  investigation  turned  to  examine  the  hypothesis  that  charge 
exchange  among  the  debris  ions  could  contribute  to  the  observed  anomalous  distribution 
of  debris.  This  investigation  considered  both  electron  capture  and  electron  loss  by  the 
debris  ions  interacting  with  the  neutral  oxygen  in  the  background.  The  hypothesis 
pursued  was  that  debris  ions  that  would  normally  be  coupled  with  the  background  plasma 
(and  therefore  not  able  to  penetrate  very  far  during  the  course  of  the  simulation)  would 
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travel  further  (radially  from  the  explosion)  with  eharge  exchange  as  a  portion  of  the  ions 
would  be  neutralized  and  therefore  not  coupled  behind  the  background  plasma.  Formed 
as  fission  fragments,  these  particles  could  later  beta  decay  and  release  high-energy 
electrons  further  from  the  center  of  the  nuclear  explosion.  Columns  2  and  3  in  Table  2  list 
the  combination  of  debris  species  examined  in  the  charge  exchange  portion  of  this 
dissertation.  Debris  Species  2  (listed  in  the  third  column  of  Table  2)  shows  the  common 
fission  fragments  used  by  LLNL  in  their  separate  simulations. 

Another  characteristic  of  the  configuration  used  in  the  simulation  to  evaluate 
charge  exchange  was  the  use  of  different  atmospheric  profiles.  The  background  plasma 
(of  flash-ionized  monatomic  oxygen)  had  two  different  configurations.  The  first 
configuration  assumed  a  constant  number  density  for  the  background  plasma  throughout 
the  radial  dimension  of  the  simulation  cells.  The  second  configuration  employed  an 
exponentially-decreasing  number  density  for  the  background  plasma  of  flash  ionized 
oxygen.  These  differences  are  shown  in  the  fourth  column  in  Table  2. 

In  both  atmospheric  configurations,  there  were  two  simulations  per  debris  species 
pairing  -  one  without  charge  exchange  and  one  with  charge  exchange.  Without  charge 
exchange,  both  species  of  the  debris  plasma  coupled  with  the  background  plasma.  With 
charge  exchange,  some  of  each  species  in  the  debris  plasma  moved  through  the  shock 
front  in  the  background  plasma.  It  is  these  debris  atoms,  which  could  later  beta  decay  at 
higher  elevations  and  introduce  high-energy  electrons,  that  appear  to  account  for  some  of 
the  observed  anomalous  debris  distributions.  Table  2  summarizes  the  simulations  carried 
out  to  evaluate  the  hypothesis  of  Charge  Exchange. 
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Simulation 

Debris 
Species  1 

Debris 
Species  2 

Atmosphere 

Charge 

Exchange 

CE-1 

U+1 

U+1 

Constant 

No 

CE-2 

U+1 

U+1 

Constant 

Yes 

CE-3 

U+1 

Xe+1 

Constant 

No 

CE-4 

u+1 

Xe+1 

Constant 

Yes 

CE-5 

u+1 

Zr+1 

Constant 

No 

CE-6 

u+1 

Zr+1 

Constant 

Yes 

CE-7 

u+1 

Pb+1 

Constant 

No 

CE-8 

u+1 

Pb+1 

Constant 

Yes 

CE-9 

u+1 

Xe+1 

Exponential 

No 

CE-10 

u+1 

Xe+1 

Exponential 

Yes 

CE-11 

u+1 

Zr+1 

Exponential 

No 

CE-1 2 

u+1 

Zr+1 

Exponential 

Yes 

CE-1 3 

u+1 

Pb+1 

Exponential 

No 

CE-1 4 

u+1 

Pb+1 

Exponential 

Yes 

Table  2.  List  of  simulations  conducted  to  examine  charge  exchange. 
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VII.  RESULTS  AND  CONCLUSIONS 


A.  OVERVIEW 

The  simulations  all  used  eylindrieal  geometry  with  the  debris  loaded  on  the  center 
axis  with  an  initial  radial  velocity  set  at  2x10^  cm/s;  this  value  was  provided  by  LLNL  as 
representative  of  the  debris  velocity.  The  undisturbed  ambient  magnetic  field  was  set  at 
0.3  Gauss  in  the  z  direction.  This  geometry  frames  the  problem  of  how  streaming  charged 
particles  moved  across  the  magnetic  field  lines. 

1.  Simulations  Examining  Collisionless  Shock  Driven  by  Multiple  Ion 
Species 

Simulations  completed  to  examine  multi-species  collisionless  shock  involved  two 
separate  debris  species.  However,  to  keep  the  energy  constant,  the  two  debris  species  had 
the  same  mass  and  used  different  charge  states  to  represent  the  difference  in  their  motion. 
The  simulations  began  with  both  debris  species  having  velocities  of  2x10^  cm/s;  using  the 
same  mass  for  the  two  debris  species  reduced  the  number  of  variables  between 
simulations  down  to  one.  The  principle  under  examination  is  the  Larmor  radius.  As 
defined  in  equation  14,  the  charge  state  of  the  ion  is  the  variable  controlled  to  achieve 
different  physical  results.  The  hypothesis  tested  was  the  effect  of  the  species  with  the 
smaller  gyroradius  (represented  in  this  case  as  the  one  with  the  higher  charge  state)  on  the 
motion  of  the  species  with  the  larger  gyroradius. 

The  simulations  in  this  portion  of  the  work  all  used  uranium  ions  with  various 
charge  states.  After  testing  the  code  in  simple  scenarios  to  ensure  the  multiple  species 
were  handled  correctly,  initial  runs  included  two  different  charge  states  of  uranium  to 
represent  separate  species.  In  these  simulations,  it  was  easy  to  observe  the  correct 
behavior  of  the  different  species  of  uranium,  and  these  runs  confirmed  ZMR’s  ability  to 
handle  multiple  species. 

Following  the  initial  simulations  where  the  results  were  intuitive,  follow-on 
simulations  involved  uranium  ions  with  non-trivial  charge  states  to  analyze  the  multi¬ 
species  effects.  In  each  of  these  simulations,  one  of  the  uranium  ions  was  at  the  +1  state 
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(providing  the  longer  gyroradius)  while  the  other  species  of  uranium  had  a  charge  state  of 
+2,  +4  or  +8.  These  charge  states  are  characterized  by  ionic  motion  with  successively 
smaller  Larmor  radi. 

There  were  two  different  configurations  of  the  uranium  ions  in  terms  of  initial 
location.  In  the  first  configuration,  all  the  debris  was  collocated  at  the  beginning  of  the 
simulation.  This  caused  the  debris  to  be  co-mingled  at  the  start  of  the  simulation  and 
stream  out  radially  at  a  constant  initial  velocity.  In  the  second  configuration,  the  debris 
was  staggered  with  the  higher  charge  state  ions  on  the  outside  and  the  lower  charge  state 
ions  on  the  inside.  This  replicates  a  debris  configuration  with  the  lighter  debris  streaming 
outward  in  front  of  the  heavier  debris. 

2.  Simulations  Completed  Examining  Charge  Exchange 

Simulations  completed  to  examine  charge  exchange  involved  two  separate  debris 
species.  Buffering  problems  with  ZMR  prevented  using  more  than  two  initial  debris 
species.  The  simulations  began  with  the  debris  species  having  a  radial  velocity  of  2x10^ 
cm/s.  The  principle  examined  in  these  simulations  involved  change  exchange  with  the 
ambient  atmosphere.  As  Figure  4  shows,  the  most  prominent  gas  at  the  elevation  of 
STARFISH  is  neutral,  atomic  oxygen  so  that  is  the  catalyst  for  all  charge  exchange 
calculations.  The  simulations  always  used  uranium  +1  as  a  debris  species  and  those 
simulations  that  used  two  debris  species  included  of  the  following  additional  ions:  xenon 
+1,  zirconium  +1  or  lead  +1.  These  debris  ions  were  chosen  as  being  representative  of  the 
likely  fission  products  produced  in  nuclear  fission. 

To  evaluate  the  effects  of  charge  exchange,  the  simulations  were  first  run  with 
charge  exchange  disabled  and  the  debris  moving  under  normal  kinematics  as  charged 
particles.  The  number  densities  of  the  debris  and  background  ions  were  adjusted  to  ensure 
the  debris  species  coupled  with  the  background  ions.  Then  the  same  simulation  was 
completed  with  charge  exchange  enabled.  These  two  simulations  were  compared  and  the 
differences  related  to  charge  exchange  and  the  movement  of  the  atoms  as  neutral  particles 
allowing  for  the  examination  of  the  role  of  charge  exchange. 
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There  were  two  different  eonfigurations  of  the  neutral,  atomic  oxygen  used  as 
catalyst  for  the  charge  exchange.  This  neutral,  atomic  oxygen  represents  the  un-ionized 
portion  of  the  atmosphere  and  not  the  ionized  portion  of  the  atmosphere  that  constitutes 
the  background  plasma.  The  neutral,  atomic  oxygen  serves  only  as  a  catalyst  for  charge 
exchange  and  does  not  interact  energetically  with  the  debris  species.  Nor  are  there  any 
calculations  of  motion  for  these  catalyst  atoms.  The  first  configuration  of  neutral,  atomic 
oxygen  assumed  a  constant  number  density  in  the  radial  distance  from  the  center  of  the 
simulation.  Because  ZMR  has  symmetry  in  the  angular  direction,  no  variable  can  have  a 
dependence  on  0,  the  angle  around  the  z  axis.  Therefore,  the  constant  atmosphere  is  a 
good  model  over  a  short  simulation  grid.  The  exponentially  decaying  atmosphere 
matches  the  data  plotted  in  Figure  4.  However,  since  ZMR  takes  the  debris  which  is 
streaming  out  isotropically  and  resets  its  motion  in  one  dimensional  radially  outward,  the 
exponential  atmosphere  discounts  the  debris  which  should  be  streaming  down  into  the 
atmosphere;  these  debris  would  hit  a  denser  atmosphere  and  generate  more  collisionless 
shock.  Using  an  exponential  atmosphere  underestimates  the  amount  of  collisionless 
shock  generated.  Since  neither  the  constant  atmosphere  or  the  exponential  atmosphere  are 
perfect  models  to  use  in  ZMR,  both  atmospheric  profiles  were  simulated  to  bound  the 
results.  Both  atmospheric  profiles  allow  for  the  examination  of  the  relevant  physics  of 
collisionless  shock  generated  by  multiple  species  and  charge  exchange.  It  should  be 
stated  that  ZMR  was  chosen  for  its  relatively  short  simulation  run  times  and  accurate 
representation  of  the  kinetic  effects  on  ions  and  these  advantageous  features  over 
compensate  for  the  slight  irregularities  with  the  atmospheric  profiles. 

B.  MULTISPECIES  COLLISIONLESS  SHOCK 

The  hypothesis  of  collisionless  shock  with  multiple  debris  species  involved  the 
lighter  debris  generating  the  shock  and  accelerating  the  heavier  debris  as  it  moved 
through  that  shock  front.  This  hypothesis  was  simulated  in  two  different  geometries.  In 
the  simulations,  instead  of  varying  the  mass,  the  charge  state  was  varied  to  adjust  the 
Larmor  radius  according  to  equation  14.  Changing  the  charge  state  ensured  the  energy 
and  number  densities  in  the  simulations  remained  constant  across  the  different  computer 
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runs.  Changing  the  charge  state  produces  the  multi-species  effect  while  reducing  the 
number  of  variables  to  be  considered. 

1,  Debris  Co-Located 

Figure  18  shows  the  expansion  of  a  multiple  species  debris  plasma  with  two 
different  charges,  one  with  charge  state  +1  and  the  second  with  charge  state  +2. 
Correspondingly,  their  gyroradi  differ  by  a  factor  of  2.  This  physics  is  similar  to  the  case 
of  two  species  both  having  a  charge  state  +1  where  one  species  is  half  the  mass  of  the 
other.  Manipulating  charge  state  is  an  easier  method  to  investigate  multiple  species  while 
keeping  the  energy  and  number  densities  constant.  In  these  simulations,  neither  debris 
species  couples  into  the  background  and  both  species  turn  with  a  normal  gyroradi.  An 
acceleration  of  the  debris  would  be  observed  as  a  jump  in  radial  velocity.  In  the  phase 
space  plot  shown  in  Figure  17,  a  jump  in  radial  velocity  would  be  observed  as  a  vertical 
shift  upwards  in  the  simulation  particles’  radial  velocity.  Figure  17  shows  both  debris 
species  (U^'  and  having  smooth  curves  to  their  motion  in  phase  space  and  therefore 
not  showing  acceleration. 
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Figure  17.  Simulation  of  two  uranium  debris  speeies  with  charge  states  +1  (green)  and 
+2(red)  expanding  into  a  background  plasma  of  0+1  (blue)  plotted  as  a 
radial  velocity  [cm/s]  versus  radial  distance  [cm].  The  smooth  profile  for 
uranium  +1  indicates  no  acceleration. 

As  the  charge  state  increases  to  +4  and  above  in  subsequent  simulations,  the 
higher  charge  state  debris  (i.e.,  the  species  with  the  smaller  Larmor  radius)  couples  into 
the  background,  but  the  debris  with  charge  +1  never  does.  Using  Hewett’s  coupling 
constant  from  equation  30  and  the  densities  and  charge  states  for  the  two  debris  species 
and  the  background  species,  the  value  of  udc  is  1.35,  greater  than  the  decoupling 
threshold  of  1.0.  Therefore,  the  expectation  is  the  debris  will  couple  with  the  background. 
Figure  18  shows  the  higher  charge  state  debris  (U"^"^)  coupled  behind  the  background 
plasma  (O^')  as  expected.  Flowever,  the  lower  charge  state  debris  (U"^*)  does  not  have  an 
increase  in  radial  velocity.  The  phase  space  plot  in  Figure  18  shows  a  smooth  profile  for 
radial  velocity  for  indicating  no  acceleration. 


69 


2E+08 


i/i 

e" 

u 


1.5E+08 


lE+08 


5E+07 


-5E+07 


-lE+08 


■i:: 


1; 

■'Ok; 
.  *'■ 


0^1 

U+1 

U+4 


I 

2E+07 


l 

4E+07 


r,  radial  distance  [cm] 


Figure  18.  Simulation  of  two  uranium  debris  species  with  charge  states  +1  (green)  and 
+4  (red)  expanding  into  a  background  plasma  of  0+1  (blue)  plotted  as  a 
radial  velocity  [cm/s]  versus  radial  distance  [cm].  The  smooth  profile  for 
uranium  +1  indicates  no  acceleration. 


As  the  charge  state  for  the  debris  species  increases,  beyond  +4,  the  observed 
trends  continue  with  the  higher  charge  state  coupling  into  the  background  and  the  debris 
at  charge  state  +1  continues  to  remain  uncoupled.  Figure  19  shows  the  results  run  at  the 
highest  charge  state  of  +8,  which  corresponds  to  the  same  Larmor  radius  of  Si^*.  Again, 
the  higher  charge  state  debris  (U^*)  coupled  behind  the  background  the  background 
plasma  (O^’).  The  phase  space  plot  in  Figure  19  shows  a  smooth  profile  for  radial 
velocity  for  indicating  no  acceleration. 
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Figure  19.  Simulation  of  two  uranium  debris  speeies  with  eharge  state  +1  (green)  and 
+8  (red)  expanding  into  a  background  plasma  of  0+1  (blue)  plotted  as  a 
radial  velocity  [cm/s]  versus  radial  distance  [cm].  The  smooth  profile  for 
uranium  +1  indicates  no  acceleration. 


With  the  two  debris  species  co-located  at  the  beginning  of  the  simulation,  the 
higher  charged  debris  species  (corresponding  to  a  smaller  gyro-radii  element  at  charge 
state  +1)  did  not  generate  a  collisionless  shock  sufficient  to  accelerate  the  uranium  +1 
debris  species. 

2,  Debris  Staggered 

To  account  for  the  geometry  and  allow  the  higher  charge  debris  to  precede  the 
lower  charge  state  debris,  the  simulations  were  rerun  with  debris  staggered.  As  discussed 
in  section  Vl.D.l,  the  higher  charge  state  debris  began  the  simulation  in  a  shell  outside  of 
the  lower  charge  state  debris.  As  the  simulation  began,  all  the  debris  pushes  into  the 
charged  background  plasma.  Figure  20  shows  that  two  debris  species  with  charge  states 
+1  and  +2  expanding  into  the  background  plasma  of  O^'.  Based  on  the  decoupling 
constant  (equation  30),  these  two  debris  species  are  decoupled  from  the  background. 
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Figure  20.  Simulation  of  two  uranium  debris  speeies  with  charge  state  +1  (red)  and  +2 
(green)  expanding  into  a  background  plasma  of  0+1  (blue)  plotted  as  a 
radial  velocity  [cm/s]  versus  radial  distance  [cm].  The  smooth  profile  for 
uranium  +1  indicates  no  acceleration. 


The  results  are  similar  to  the  co-located  results  (Figure  17).  In  the  staggered 
simulation,  the  begins  outside  the  and  spends  more  time  interacting  with  the 
ambient  magnetic  field.  The  debris  follows  behind  the  and  passes  through  the 
collisionless  shock  generated  by  the  However,  in  Figure  20,  the  profile  of  the  radial 
velocity  for  remains  smooth  showing  no  acceleration. 

As  the  charge  state  of  the  outside  debris  species  increases  to  +4,  the  results  mirror 
the  co-located  scenarios.  Figure  21  shows  the  coupled  inside  the  background  plasma 
with  the  not  trapped  inside  the  background. 
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Figure  21.  Simulation  of  two  uranium  debris  species  with  charge  state  +1  (red)  and  +4 
(green)  expanding  into  a  background  plasma  of  0+1  (blue)  plotted  as  a 
radial  velocity  [cm/s]  versus  radial  distance  [cm].  The  smooth  profile  for 
uranium  +1  indicates  no  acceleration. 


As  shown  in  Figure  21,  has  a  smooth  profile  in  radial  velocity  indicating  no 
acceleration  from  the  collisionless  shock  formed  by  the  Similarly,  as  the  charge  state 
of  the  debris  species  increases  to  +8,  the  inside  debris  species  remains  decoupled  as  seen 
in  Figure  22.  Again,  the  U^*,  which  represents  a  larger  gyroradius,  maintains  a  smooth 
profile  for  the  radial  velocity  in  phase  space.  Figure  22  shows  no  indication  of 
acceleration  due  to  the  collisionless  shock  formed  by  the  (representing  a  smaller 
gyroradius). 
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Figure  22.  Simulation  of  two  uranium  debris  speeies  with  charge  state  +1  (red)  and  +8 
(green)  expanding  into  a  background  plasma  of  0+1  (blue)  plotted  as  a 
radial  velocity  [cm/s]  versus  radial  distance  [cm].  The  smooth  profile  for 
uranium  +1  indicates  no  acceleration. 


This  simulation  did  not  support  the  hypothesis  that  the  smaller  gyroradius  debris 
would  drive  the  collisionless  shock  and  accelerate  the  larger  gyroradius  debris.  The 
simulation  suggests  the  debris  species  interact  independently  with  the  ambient  magnetic 
field  and  the  background  plasma.  Even  with  the  debris  staggered  at  the  beginning  of  the 
simulation,  the  smaller  gyroradius  debris  did  not  form  a  collisonless  shock  large  enough 
to  accelerate  the  larger  gyroradii  debris. 

C.  CHARGE  EXCHANGE 

The  effects  from  charge  exchange  were  examined  using  the  same  initial 
conditions  for  a  multi-species  simulation  where  the  debris  coupled  with  the  shock  front 
and  remained  inside  the  background  ions.  That  is,  as  the  simulation  progressed  all  the 
debris  stayed  behind  the  shock  front  and  expended  all  its  energy  driving  the  shock.  These 
same  initial  conditions  were  simulated  with  charge  exchange  enabled.  Since  the  initial 
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simulation  had  100%  of  the  debris  eoupled  inside  the  baekground  plasma,  any  debris  now 
ahead  of  the  shoek  front  is  attributed  to  the  effeets  of  eharge  exehange. 

1.  Sensitivity 

The  sensitivity  of  the  number  of  neutral  atoms  produeed  to  the  eross  seetions  was 
examined.  The  sensitivity  of  the  eross  seetion  for  both  eleetron  eapture  and  eleetron  loss 
was  examined  independently.  First  the  eross  seetion  for  eleetron  eapture  was  sealed  by  a 
faetor  of  two,  then  five  and  finally  ten.  The  total  number  of  neutral  atoms  at  the  end  of 
the  simulation  sealed  up  aeeordingly.  Seeondly,  the  eross  seetion  for  eleetron  loss  was 
sealed  by  a  faetor  of  two,  then  five  and  finally  ten.  Then  the  total  number  of  neutral 
atoms  at  the  end  of  the  simulation  sealed  down  aeeordingly.  As  expeeted,  the  total 
number  of  neutrals  produeed  direetly  seales  with  the  values  of  the  eross  seetions  for  both 
eleetron  eapture  and  eleetron  loss  eomputed  by  the  semi-empirieal  formula. 

The  sensitivity  of  the  eross  seetions  is  important  beeause  the  semi-empirieal 
formula  does  not  exaetly  mateh  the  experimental  data.  If  better  semi-empirieal  formulas 
for  the  eross  seetions  of  the  elements  of  interest  are  found,  they  eould  easily  be 
implemented  in  the  simulation.  The  values  determined  from  the  semi-empirieal  here  are 
used  to  demonstrate  the  relevant  physies  of  eharge  exehange. 

2,  Effects  of  Collisionless  Shock  on  Charge  Exchange 

The  impaet  of  eollisionless  shook  on  eharge  exehange  was  examined.  The 
simulation  was  run  without  a  oharged  baekground  plasma  using  neutral  hydrogen  as  the 
baekground  material.  The  neutral  oxygen  profile  was  also  loaded  and  eharge  exehange 
enabled.  This  oaused  the  expanding  debris  of  to  stream  radial  outward  and  without 
the  oharged  baekground  plasma  there  is  no  ooupling  between  the  debris  plasma  and  the 
baekground  so  there  is  no  kinetio  energy  transferred  to  the  baekground  plasma  by  shook 
formation.  As  the  oharged  baekground  passes  through  the  neutral  oxygen  profile,  eharge 
exehange  ooours  but  without  any  impaet  from  eollisionless  shook  (whioh  is  absent  due  to 
the  laok  of  a  eharge  baekground).  The  total  number  of  neutral  uranium  atoms  as  a 
funotion  of  time  is  shown  in  Figure  23.  The  shape  of  the  results  plotted  in  Figure  23 

mateh  the  expeotation  from  equation  49  and  theoretioal  oaloulation  plotted  in  Figure  15. 
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The  importance  of  Figure  23  is  the  value  of  the  total  number  of  real,  neutral  uranium 
atoms  produced  in  the  simulation. 


Without  Collisionless  Shook 


Figure  23.  Simulation  of  uranium+1  expanding  into  a  background  of  neutral  hydrogen 
without  generating  a  collisionless  shock  and  with  charge  exchange  enabled. 

The  same  conditions  were  used  in  another  simulation  however  in  this  case  the 
background  plasma  was  switched  from  neutral  hydrogen  back  to  O^'.  The  charged 
background  and  the  ambient  magnetic  field  interact  with  the  expanding  debris  ions  (U^*) 
generating  a  collisionless  shock.  Again,  the  neutral  oxygen  profile  remained  on  the  grid 
and  charge  exchange  was  enabled.  The  difference  between  the  two  simulations  was  only 
the  formation  of  the  collisionless  shock.  The  simulation  results  of  the  number  of  real, 
neutral  uranium  atoms  plotted  as  a  function  of  time  is  shown  in  Figure  24.  Again,  the  data 
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plotted  in  Figure  24  matehes  the  form  expeeted  from  equation  49  and  the  late  time  value 
of  the  total  number  of  neutral  uranium  atoms  produeed  is  important. 


Figure  24.  Simulation  of  uranium+1  expanding  into  a  baekground  of  oxygen+1  while 
generating  a  collisionless  shock  and  with  charge  exchange  enabled.  Plot  of 
number  of  real,  neutral  uranium  atoms  as  a  function  of  time. 

The  simulations  were  similar  except  for  the  collisionless  shock.  Both  began  with 
the  same  energy  and  only  the  collisionless  shock  moved  energy  from  the  debris  to  the 
background  plasma.  Both  simulations  reached  a  steady  state  of  neutral  uranium  in  about 
0.06  seconds  however,  the  steady  state  values  for  the  two  simulations  were  different.  The 
simulation  without  collisionless  shock  reached  a  steady  state  value  of  1.26x10^^  atoms  of 
neutral  uranium  and  the  simulation  with  collisionless  shock  reached  a  steady  state  value 
of  8x10^^  atoms  of  neutral  uranium.  While  the  absence  of  a  collisionless  shock  caused  an 
increase  of  33%  in  the  number  of  neutrals,  the  net  change  in  velocity  with  and  without 

collisionless  shock  was  only  1%.  Without  collisionless  shock  and  no  mechanism  for 
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energy  exehange  between  the  debris  ions  and  the  baekground,  the  debris  ions  lost  no 
energy  or  veloeity  to  the  baekground  and  maintained  their  veloeity  of  2x10^  ern/s.  With 
eollisionless  shook  and  a  meohanism  to  exehange  energy  between  the  debris  ions  and  the 
baekground  ions,  the  average  speed  of  the  neutral  partioles  was  1.98x10*  om/s.  The 
probability  of  oharge  exehange  is  oaloulated  using  equation  34  and  inoludes  oross  seotion, 
speed,  size  of  the  time  step,  and  number  density  of  the  target  (neutral  oxygen). 
Additionally,  the  oross  seotion  is  energy  dependent  and  is  influenoed  by  the  debris  ion’s 
speed.  A  deorease  of  1%  in  veloeity  ohanges  the  oross  seotion  for  eleotron  oapture  by  5% 
and  ohanges  the  oross  seotion  for  eleotron  loss  by  3%.  Both  with  and  without 
eollisionless  shook,  the  size  of  the  time  step  and  the  number  density  profile  of  the  target 
atoms  is  the  same.  Neither  the  size  of  the  time  step  nor  the  number  density  profile  of  the 
target  atoms  aooounts  for  the  differenoe  between  the  numbers  of  neutrals  produoed  for  the 
simulations  with  and  without  eollisionless  shook. 

Collisionless  shook  introduoes  an  eleotrio  field  in  the  angular  direotion  that  turns 
the  debris  ions  from  their  original  radial  path.  The  angular  eleotrio  field  is  small  and  the 
overall  effeot  on  the  speed  is  small;  henoe  the  only  1%  differenoe  in  the  speeds.  However, 
the  speed  is  the  veotor  sum  of  the  oomponents  of  the  veloeity  veotor.  While  the 
magnitude  of  the  veloeity  is  the  same,  the  veloeity  oomponents  have  shifted  with  a  larger 
theta  veloeity.  Due  to  the  turning,  the  debris  ions  are  no  longer  streaming  radially  through 
an  exponentially  deoreasing  atmosphere  with  a  deoreasing  number  density  of  target  atoms 
(neutral  oxygen).  Instead,  as  the  debris  ions  turn  in  the  presenoe  of  a  transverse  eleotrio 
field,  the  debris  ions  move  through  a  higher  number  density  of  target  atoms  and  undergo 
more  oharge  exehange  before  the  oharge  exehange  “freezes  out.” 

3.  Constant  and  Exponential  Atmosphere 

As  desoribed  in  Chapter  V.C,  the  simulations  were  run  with  both  a  oonstant  and 
exponential  atmosphere  and  there  are  different  results  in  slip  through  for  eaoh 
atmosphere.  The  different  atmospherio  profiles  produoed  dramatioally  different  results  in 
the  amount  of  slip  through.  In  oylindrioal  geometry,  neither  atmospherio  representation  is 
a  oorreot  representation  of  the  aotual  atmosphere;  however,  this  work  shows  there  is  a 
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clear  difference  in  the  results  and  the  appropriate  atmospheric  model  should  be  used. 
Future  work  in  three  dimensions  should  include  more  accurate  atmospheric  model. 

4,  Slip  through  with  Constant  Atmosphere 

The  results  described  in  this  section  are  from  simulations  using  a  constant 
atmosphere  for  both  the  background  plasma  (O^*)  and  the  neutral  oxygen  background 
(0"«). 


a.  Uranium  Only 

These  simulations  show  the  effects  of  charge  exchange,  demonstrate  the  results 
and  describe  the  process  that  will  be  repeated  when  a  second  debris  species  is  added  in 
later  simulations.  The  simulation  included  one  Uranium  debris  species  with  a  charge  state 
of  +1.  In  the  first  simulation  there  was  no  charge  exchange.  The  debris  began  with  an 
initial  radial  velocity  of  2x10*  cm/s  which  corresponds  to  a  kinetic  energy  representative 
of  0.29  MT  yield.  Without  charge  exchange,  the  debris  couples  with  the  background  ions 
as  seen  in  Figure  25. 


Figure  25.  Simulation  of  Uranium+1  (red)  expanding  into  background  plasma  of 

Oxygen+1  (blue)  in  a  constant  atmosphere  without  charge  exchange. 
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The  simulation  was  rerun  again  with  only  one  debris  species  of  but  with 
charge  exchange  enabled.  During  this  simulation,  particles  experienced  charge  exchange 
of  both  electron  capture  and  electron  loss.  All  other  inputs  were  kept  the  same.  Figure  26 
shows  different  results  from  charge  exchange  with  some  of  the  uranium  atoms  slipping 
ahead  of  the  shock  front  in  the  background  plasma.  The  comparison  between  Figures  25 
and  26  show  the  clear  difference  the  effects  from  simulations  run  with  charge  exchange, 
specifically  how  much  further  down  the  simulation  grid  the  debris  travels  and  how  some 
debris  clearly  is  decoupled  from  the  background  plasma. 
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Figure  26. 


Simulation  of  Uranium+1  expanding  into  background  plasma  of  Oxygen+1 
in  a  constant  atmosphere  with  charge  exchange. 


Enabling  charge  exchange  allowed  some  debris  to  slip  ahead  of  the  shock  front. 
Since  the  initial  simulation  without  charge  exchange  resulted  in  all  the  debris  behind  the 
shock  front,  this  makes  an  easy  comparison  to  see  the  effects  of  charge  exchange.  Charge 
exchanged  allowed  1.02%  of  the  debris  to  slip  ahead  of  the  shock  front  with  most  of  the 
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slip-through  debris  being  in  the  neutral  charge  state.  A  comparison  of  the  results  between 
the  two  simulations  with  as  the  debris  species  is  shown  in  Table  3  with  the  last  two 
columns  showing  the  amount  of  real  atoms  that  slipped  ahead  of  the  magnetic  shock  front 
at  the  completion  of  the  simulation. 


Total#  of 
ions  (t=0) 

Total#  of 
atoms  (t=0) 

Total#  of 
ions  ahead  of 
shock  front 
(t=5000) 

Total#  of 
ions  ahead  of 
shock  front 
(t=5000) 

Without  Charge 
Exchange 

1.70x10^^ 

0 

0 

N/A 

With  Charge 
Exchange 

1.70x10^^ 

0 

5.18x10^* 

1.69x10^^ 

Table  3.  Comparison  of  slip  through  data  without  and  with  charge  exchange  for 

uranium. 


Figures  25  and  26  illustrate  significant  differences  in  simulation  outcomes 
(though  it  should  be  noted  that  Figure  26  includes  more  points  plotted  than  Figure  25). 
Both  simulations  have  the  same  number  of  particles  and  energy.  However,  the  use  of  the 
charge  exchange  function  in  ZMR  creates  new  simulation  particles  and  reduces  the  real- 
particle-per-simulation  particle  value.  Therefore,  Figure  26  has  more  points  plotted  but 
each  point  represents  fewer  real  particles  in  the  simulation.  The  same  effect  will  occur  for 
all  simulations  involving  charge  exchange  and  all  those  graphs  will  have  more  points 
plotted  at  the  completion  of  the  simulation  than  the  corresponding  runs  without  charge 
exchange. 

b.  Uranium  and  Xenon 

The  next  set  of  simulations  includes  uranium  and  xenon  both  at  charge  state  +1. 
Again,  the  debris  began  the  simulation  with  an  initial  velocity  of  2x10*  cm/s  in  the  radial 
direction.  The  kinetic  energy  representative  of  a  0.29  MT  yield.  At  the  completion  of  the 
simulation  with  charge  exchange  disabled,  the  debris  is  coupled  to  the  background 
plasma.  Figure  27  shows  the  debris  of  U^'  (light  blue)  and  Xe"^'  (red)  behind  the 
background  plasma  O^'  (blue).  Without  charge  exchange,  100%  of  the  debris  couples  to 
the  background  plasma  and  none  of  the  debris  gets  ahead  of  the  shock  front. 
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r,  radial  distance  [cm] 

Figure  27.  Simulation  involving  Uranium+1  (light  blue)  and  Xenon+1  (red)  expanding 
into  background  plasma  of  Oxygen+1  with  a  constant  atmosphere  without 

charge  exchange. 

The  simulation  was  run  again  with  charge  exchange  enabled.  All  the  other 
parameters  were  kept  constant.  The  simulation  allowed  the  ions  to  neutralize  to 
and  later  re-ionize  back  to  U^*.  The  same  was  true  to  the  Xe^*  ions.  Figure  28  shows  the 
results  of  the  simulation  at  its  completion  with  some  of  the  debris  ahead  of  the  shock 
front.  The  comparison  between  Figures  27  and  28  show  the  difference  charge  exchange 
contributes  to  the  results  with  two  species. 
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Figure  28.  Simulation  involving  Uranium+1  (light  blue)  and  Xenon+1  (gold) 
expanding  into  a  baekground  plasma  of  Oxygen+1  with  a  eonstant 
atmosphere  and  charge  exchange  enabled. 


Without  charge  exchange,  all  the  debris  coupled  inside  the  background  plasma 
and  none  was  in  front  of  the  shock  front.  With  charge  exchange,  some  of  the  debris  was 
able  to  move  in  front  of  the  shock  front  and  beyond  the  background  plasma.  Charge 
exchange  resulted  in  0.33%  of  the  debris  ending  in  front  of  the  shock  front.  Table  4 
shows  the  results  of  the  two  simulations  containing  U  and  Xe  as  debris  with  charge 
exchange  disabled  and  then  again  with  charge  exchange  enabled.  The  last  four  columns 
of  Table  4  show  the  number  of  atoms  that  moved  beyond  the  shock  front  at  the 
completion  of  each  simulation. 
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Total  # 
of 
ions 
(t=0) 

Total  # 
of 
ions 
(t=0) 

Total  # 
of  U^® 
atoms 
(t=0) 

Total  # 
of  Xe^® 
atoms 
(t=0) 

Total  # 

of 

ions 

ahead 

of  shock 

front 

(t=5000) 

Total  # 

of  Xe^‘ 

ions 

ahead 

of  shock 

front 

(t=5000) 

Total  # 

of  U^® 

ions 

ahead 

of  shock 

front 

(t=5000) 

Total  # 

of  Xe^® 

ions 

ahead 

of  shock 

front 

(t=5000) 

Without 

Charge 

Exchange 

1.41x10^^ 

2.81x10^'’ 

0 

0 

0 

0 

NA 

NA 

With 

Charge 

Exchange 

1.41x10^^ 

2.81x10“ 

0 

0 

1.25x10“ 

7.94x10“ 

4.34x10“ 

6.49x10“ 

Table  4.  Comparison  of  slip  through  data  for  uranium  and  xenon  without  and  with 

eharge  exehange  in  a  eonstant  atmosphere. 


c.  Uranium  and  Zirconium 

The  next  set  of  simulations  ineluded  both  uranium  and  zireonium  as  debris 
speeies  with  eharge  state  +1.  Again,  the  debris  began  the  simulation  with  a  veloeity  of 
2x10*  em/s  in  the  radial  direetion.  The  kinetie  energy  of  the  debris  was  representative  of  a 
0.29  MT  yield.  At  the  eompletion  of  the  simulation  with  eharge  exehange  disabled,  the 
debris  had  eoupled  inside  the  background.  Figure  29  shows  the  debris  species  coupled  to 
the  background  plasma  at  the  completion  of  the  simulation. 
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Figure  29.  Simulation  involving  Uranium+1  (green)  and  Zirconium+1  (red)  expanding 
into  a  background  plasma  of  0+1  (blue)  with  a  constant  atmosphere  without 

charge  exchange. 

The  simulation  was  rerun  with  all  the  same  parameters  and  charge  exchange 
enabled.  With  charge  exchange,  some  of  the  debris  was  neutralized  and  slipped  ahead  of 
the  shock  front.  Figure  30  shows  the  results  of  the  simulation  with  charge  exchange 
enabled.  A  comparison  of  Figures  29  and  30  show  the  difference  charge  exchange 
contributes  to  the  results  with  two  species. 


85 


[cm/s] 


2E+08 


lE+08 


-lE+08 


0 


u 

Zr 

zr° 

0*1 


i_I 

2E+07  4E+07 


r,  radial  distance  [cm] 


i 

6E+07 


Figure  30.  Simulation  involving  Uranium+1  (light  blue)  and  Zirconium+1  (gold) 
expanding  into  a  baekground  plasma  of  0+1  (blue)  with  a  constant 
atmosphere  and  charge  exchange. 


Enabling  charge  exchange  in  the  simulation  allowed  1.1%  of  the  debris  to  slip 
ahead  of  the  shock  front.  Table  5  shows  the  results  of  both  simulations  involving  U  and 
Zr  at  the  completion  of  the  simulations  with  charged  exchange  disabled  and  then  again 
with  charge  exchange  enabled.  The  last  four  columns  of  Table  5  show  the  difference  in 
the  number  of  atoms  that  slip  beyond  the  shock  front  at  the  completion  of  the  simulations. 
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Total  # 
of 
ions 
(t=0) 

Total  # 
of 
ions 
(t=0) 

Total 

#  of 

^+0 

atoms 

(t=0) 

Total 
#  of 

zC” 

atoms 

(t=0) 

Total  # 

of 

ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 
of  Zr^‘ 
ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 
of  U^® 
ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 
of  Zr^® 
ions 

ahead  of 
shock 
front 
(t=5000) 

Without 

Charge 

Exchange 

1.41x10^^ 

3.7x10“ 

0 

0 

0 

0 

NA 

NA 

With 

Charge 

Exchange 

1.41x10^’ 

3.7x10“ 

0 

0 

5.74x10“ 

3.31x10“ 

1.40x10“ 

2.31x10“ 

Table  5.  Comparison  of  slip  through  data  without  and  with  charge  exchange  for 

uranium  and  zirconium  in  a  constant  atmosphere. 


d.  Uranium  and  Lead 

The  next  set  of  simulations  includes  both  uranium  and  lead  both  at  charge  state 
+1.  Again,  the  debris  began  the  simulation  with  a  radial  velocity  of  2x10*  cm/s  and  the 
total  kinetic  energy  of  the  debris  was  representative  of  a  0.29  MT  yield.  At  the 
completion  of  the  simulation  with  charge  exchange  disabled,  all  the  debris  had  coupled 
behind  the  shock  front.  Figure  31  shows  the  results  at  the  completion  of  the  simulation 
with  charge  exchange  disabled. 
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Figure  31.  Simulation  involving  Uranium+1  (green)  and  Lead+1  (red)  expanding  into  a 
background  of  Oxygen+1  (blue)  with  a  constant  atmosphere  without  charge 

exchange. 

The  simulation  was  rerun  with  all  the  same  parameters  with  charge  exchanged 
enabled.  Figure  32  below  shows  the  results  of  the  simulation  with  charge  exchange 
enabled.  At  the  completion  of  the  simulation,  some  of  the  debris  was  in  front  of  the  shock 
front.  The  comparison  between  Figures  31  and  32  shows  the  impact  of  charge  exchange 
on  the  simulation  results. 
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Figure  32.  Simulation  with  Uranium  +1  (light  blue)  and  Lead+1  (gold)  expanding  into 
a  background  plasma  of  Oxygen+1  (blue)  in  a  constant  atmosphere  with 

charge  exchange. 


Charge  exchange  enabled  1%  of  the  debris  to  slip  ahead  of  the  shock  front.  Table 
6  shows  a  comparison  of  the  results  of  uranium  and  lead  with  charge  exchange  disabled 
and  then  with  charge  exchange  enabled.  The  last  four  columns  show  the  total  number  of 
atoms  ahead  of  the  shock  front  at  the  completion  of  the  simulation. 
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Total  # 
of 
ions 
(t=0) 

Total  # 
of 
ions 
(t=0) 

Total 

#  of 

^+0 

atoms 

(t=0) 

Total 
#  of 
Pb^® 
atoms 
(t=0) 

Total  # 

of 

ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 
of  Pb^‘ 
ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 
of  U^® 
ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 
of  Pb^® 
ions 

ahead  of 
shock 
front 
(t=5000) 

Without 

Charge 

Exchange 

1.41x10^^ 

1.8x10^'’ 

0 

0 

0 

0 

NA 

NA 

With 

Charge 

Exchange 

1.41x10^’ 

1.8x10“ 

0 

0 

6.31x10“ 

4.41x10“ 

1.44x10“ 

1.31x10“ 

Table  6.  Comparison  of  slip  through  data  without  and  with  charge  exchange  for 

uranium  and  lead  in  a  constant  atmosphere. 


5.  Slip  through  with  Exponential  Atmosphere 

a.  Uranium  Only 

The  initial  simulation  examining  charge  exchange  with  the  assumption  of  an 
exponential  atmosphere  involved  only  Uranium  at  charge  state  +1.  The  initial  debris 
velocity  was  set  at  2x10*  cm/s  and  the  number  density  of  was  set  at  7x10^  part/cm*. 
This  was  set  to  have  the  single  species  couple  inside  the  background  plasma.  It  also 
corresponds  to  kinetic  energy  of  a  0.25  MT  yield.  These  initial  conditions  ensure  the 
debris  plasma  will  couple  to  the  background  plasma  to  show  the  impact  of  charge 
exchange  on  the  final  spatial  distribution  of  the  debris  ions.  It  is  motion  of  these  debris 
ions,  which  will  later  beta  decay  and  inject  high-energy  electrons,  that  is  the  foundation 
of  this  research.  Figure  33  shows  the  results  at  the  completion  of  the  simulation  with 
charge  exchange  disabled.  Here  all  the  debris  is  coupled  with  the  background  plasma. 
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Figure  33.  Simulation  involving  Uranium  +1  (Red)  into  Oxygen  +1  (Blue)  without 
eharge  exchange  where  the  debris  couples  into  background. 

The  same  simulation  was  run  with  charge  exchange  enabled.  All  the  other 
parameters  were  kept  the  same.  All  the  debris  ions  initially  start  at  charge  state  +1.  Figure 
34  shows  the  phase  space  distribution  at  the  end  of  the  simulation.  The  debris  ions  of 
are  plotted  in  green  and  the  neutral  uranium  atoms  (U^*’)  are  plotted  in  red.  The 
comparison  between  Figures  33  and  34  illustrates  the  impact  of  charge  exchange  on  the 
results  of  the  simulation.  Figure  34  shows  that  some  of  the  debris  has  slipped  through  the 
background  plasma. 
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Figure  34.  Simulation  involving  Uranium  +1  (Green)  into  Oxygen  +1  (Blue)  with 
eharge  exehange  forming  neutral  Uranium  (Red)  where  some  debris  slips 

ahead  of  shook  front. 


With  eharge  exehange  enabled,  6%  of  the  debris  speoies  (uranium)  gets  ahead  of 
the  shook  front.  The  magnetio  field  shook  front  in  the  z  oomponent  was  defined  as  0.4 
Gauss.  The  total  number  of  real  partioles  in  front  of  the  shook  front  was  oaloulated  by 
ZMR.  Without  eharge  exehange,  there  is  no  meohanism  to  produoe  the  U^°  atoms  so  they 
oould  not  be  in  front  of  the  shook  front  in  the  first  simulation  (seoond  row  of  table).  Table 
7  shows  the  oomparison  between  the  simulations  with  eharge  exehange  disabled  and  with 
eharge  exehange  enabled.  The  last  two  oolumns  show  the  number  of  atoms  that  slip 
ahead  of  the  shook  front. 
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Total  #  of 
ions  (t=0) 

Total  #  of 
atoms  (t=0) 

Total  #  of 
ions  ahead  of 
shock  front 

(t=5000) 

Total  #  of 
ions  ahead  of 
shock  front 

(t=5000) 

Without  Charge 
Exchange 

1.37x10^’ 

0 

0 

NA 

With  Charge 
Exchange 

1.37x10^’ 

0 

8.35x10^^ 

8.78x10^^ 

Table  7.  Comparison  of  slip  through  data  without  and  with  charge  exchange  for 

uranium. 


b.  Uranium  and  Xenon 

o 

Again,  all  the  debris  began  the  simulation  with  a  speed  of  2x10  [cm/s]  and  the 
number  densities  were  adjusted  so  the  system  had  kinetic  energy  representative  of  a  0.25 
MT  yield.  The  simulation  was  first  completed  without  charge  exchange  and  the  debris 
clearly  coupled  behind  the  shock  front  at  the  end  of  the  simulation.  Figure  35  shows  both 
of  the  debris  species  coupled  behind  the  background  plasma  at  the  completion  of  the 
simulation. 


Figure  35.  Simulation  involving  Uranium  +1  (green)  and  Xenon  +1  (red)  into  Oxygen 
+1  (blue)  without  charge  exchange  where  debris  couples  into  background. 
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The  simulation  was  run  again  with  all  the  same  parameters  and  charge  exchange 
was  enabled.  Figure  36  shows  the  phase  space  distribution  of  the  debris  at  the  completion 
of  the  simulation.  The  comparison  between  Figures  35  and  36  show  the  impact  of  charge 
exchange  in  the  simulation. 


Figure  36.  Simulation  involving  Uranium+1  (light  blue)  and  Xenon+1  (red)  into 
Oxygen+1  (blue)  with  charge  exchange  forming  neutral  Uranium  (green) 
and  neutral  Xenon  (gold)  where  some  debris  slips  through  the  background. 

ZMR  calculated  the  portion  of  the  debris  ahead  of  the  shock  front  which  was  set 
at  Bz=0.4  Gauss.  With  charge  exchange  enabled,  10.7%  of  the  debris  slipped  through  the 
shock  in  the  background  plasma.  The  results  are  shown  in  Table  8  below.  The  last  four 
columns  show  the  number  of  atoms  ahead  of  the  shock  front  at  the  completion  of  the 
simulation.  Without  charge  exchange,  there  is  no  mechanism  for  the  simulation  to 
produce  or  Xe^**  so  they  could  not  be  in  front  of  the  shock  front  at  the  5000  time  step 
in  the  first  simulation  (second  row  of  table). 
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Total# 

ofU-> 

ions 

(t=0) 

Total# 

ofXe*> 

ions 

(t=0) 

Total# 

ofU-« 

atoms 

(t=0) 

Total# 

ofXe-« 

atoms 

(t=0) 

Total# 

off* 

ions 

ahead of 
shock 
front 
(t=5000) 

Total# 

ofXe-» 

ions 

ahead of 
shock 
front 
(t=5000) 

Total# 

ofL7-« 

ions 

aheadof 

shock 

front 

(t=5000) 

Total# 

ofXe-« 

ions 

aheadof 

shock 

front 

(t=5000) 

tMthout 

Charge 

Exchange 

1.24xl(F 

2.47xl(P* 

0 

0 

0 

0 

NA 

NA 

With 

Charge 

Exchange 

1.24xl(F 

2.47x10* 

0 

0 

1.84x10* 

7,92x10* 

8.05x10“ 

1.43x10* 

Table  8.  Comparison  of  slip  through  data  without  and  with  charge  exchange  for 

uranium  and  xenon. 


c.  Uranium  and  Zirconium 

o 

Again,  the  debris  began  the  simulation  with  a  speed  of  2x10  [cm/s]  and  the 
number  densities  were  adjusted  so  the  system  had  kinetic  energy  representative  of  0.25 
MT  of  yield.  The  simulation  was  initially  run  without  charge  exchange  enabled.  Figure 
37  shows  the  debris  clearly  coupled  behind  the  shock  front  at  the  end  of  the  simulation. 
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[cm/s] 


r,  radial  distance  [cm] 


Figure  37.  Simulation  involving  Uranium  +1  (green)  and  Zirconium  +1  (red)  into 
Oxygen  +1  (blue)  without  charge  exchange  where  debris  couples  into 

background. 

The  simulation  was  run  again  with  all  the  same  parameters  and  charge  exchange 
was  enabled.  Figure  38  shows  the  phase  space  disposition  of  the  debris  at  the  completion 
of  the  simulation  with  charge  exchange.  The  comparison  between  Figures  37  and  38 
shows  the  impact  on  the  simulation  results  with  charge  exchange  enabled. 
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Figure  38.  Simulation  involving  Uranium  +1  (Blue)  and  Zireonium  +1  (Yellow)  into 
Oxygen  +1  (Blue)  with  eharge  exchange  where  some  debris  slips  through 

the  background. 

At  the  completion  of  the  simulation  with  charge  exchange  enabled,  2.3%  of  the 
total  debris  had  slipped  through  the  shock  front,  and  nearly  all  that  debris  was  neutral. 
There  were  at  least  six  orders  of  magnitude  between  the  number  of  neutral  and  charged 
constituents  in  the  debris  ahead  of  the  shock  front.  The  results  are  shown  in  Table  9 
below.  The  last  four  columns  show  the  number  of  atoms  ahead  of  the  shock  front  at  the 
completion  of  the  simulation.  Without  charge  exchange,  there  is  no  mechanism  for  the 
simulation  to  produce  or  Zr^°  so  they  could  not  be  in  front  of  the  shock  front  at  the 
5000  time  step  in  the  first  simulation  (second  row  of  table). 
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Total  # 
ofU*i 

ions 

(t=0) 

Total  # 
ofZi-‘ 
ions 
(t=0) 

Total  # 
ofU^ 

atoms 

(t=0) 

Total  # 
ofZr-« 

atoms 

(t=0) 

Total  # 
ofU*i 

ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 
ofZr-i 

ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 
ofU^o 

ions 

ahead  of 
shock 
front 
(t=5000) 

Total  U 
ofZr-« 
ions 

ahead  of 
shock 
front 
(t=5000) 

Without 

Charge 

Exchange 

1.24x10^^ 

3.55xl(P‘ 

0 

0 

0 

0 

NA 

NA 

With 

Charge 

Exchange 

1.24x10^7 

3.55xl(F‘ 

0 

0 

1.51x1015 

9.76x1011 

4.31x101^1 

6.79x101^1 

Table  9.  Comparison  of  slip  through  data  without  and  with  charge  exchange  for 

uranium  and  xenon. 


d.  Uranium  and  Lead 

Q 

Again,  the  debris  began  the  simulation  with  a  speed  of  2x10  [cm/s]  and  the 
number  densities  were  adjusted  so  the  system  had  kinetic  energy  representative  of  0.25 
MT  of  yield.  The  simulation  was  initially  run  without  charge  exchange.  Figure  39  shows 
the  debris  clearly  coupled  behind  the  shock  front  at  the  end  of  the  simulation. 
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[cm/s] 


r,  radial  distance  [cm] 

Figure  39.  Simulation  involving  Uranium  +1  (Green)  and  Lead  +1  (red)  into  Oxygen 
+1  (Blue)  without  eharge  exehange  where  debris  eouples  into  baekground. 

The  simulation  was  run  again  with  all  the  same  parameters  and  charge  exchange 
was  enabled.  Figure  40  shows  the  phase  space  disposition  of  the  debris  at  the  completion 
of  the  simulation  with  charge  exchange.  The  comparison  between  Figures  39  and  40 
shows  the  impact  on  the  simulation  results  with  charge  exchange  enabled. 
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r,  radial  distance  [cm] 


Figure  40.  Simulation  involving  Uranium  +1  (Blue)  and  Lead  +1  (Yellow)  into 
Oxygen  +1  (Blue)  with  charge  exchange  where  debris  slips  through 

background. 


At  the  completion  of  the  simulation,  6.4%  of  the  total  debris  had  passed  through 
the  shock  front,  and  nearly  all  that  debris  was  neutral.  There  was  a  two  order  of 
magnitude  difference  between  the  neutral  and  charged  constituents  of  the  debris  ahead  of 
the  shock  front.  The  results  are  shown  in  Table  10  below.  The  last  four  columns  show  the 
number  of  atoms  ahead  of  the  shock  front  at  the  completion  of  the  simulation.  Without 
charge  exchange,  there  is  no  mechanism  for  the  simulation  to  produce  or  so  they 
could  not  be  in  front  of  the  shock  front  at  the  5000  time  step  in  the  first  simulation 
(second  row  of  table). 
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Total  # 
oft*i 

ions 

(t=0) 

Total  # 
ofPb*i 

ions 

(t=0) 

Total  # 
ofU*« 

atoms 

(t=0) 

Total  # 
of  Pb*« 

atoms 

(t=0) 

Total  # 

ofU*‘ 

ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 
ofPb^i 

ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 

ofU*« 

ions 

ahead  of 
shock 
front 
(t=5000) 

Total  # 
of  Pb*« 
ions 

ahead  of 
shock 
front 
(t=5000) 

Without 

Charge 

Exchange 

1.24xl(P' 

1.55xl0^« 

0 

0 

0 

0 

NA 

NA 

With 

Charge 

Exchange 

1.24xl(F' 

1.55xl0^« 

0 

0 

1.89x10=3 

1.39x10=3 

7.92x10=3 

9.86x10='' 

Table  10.  Comparison  of  slip  through  data  without  and  with  charge  exchange  for 

uranium  and  lead. 


6,  Calculations  of  Neutral  Yield  to  Match  Measured  Flux 

Assuming  the  fission  fraction  of  the  total  yield  is  0.5  and  each  fission  produces 

26 

two  fission  fragments,  STARFISH  would  produce  2x10  fission  fragments.  Assuming 
four  beta  decays  per  fission  fragment,  this  corresponds  to  8x10  electrons  from  beta 
decay.  Hess  states  the  total  number  of  trapped  electrons  between  1.75  Earth  radii  (Re)  and 
2.25  Re  is  1.2x10^^  electrons.  Hess’  number  of  trapped  electrons  corresponds  to  15%  of 
the  total  electrons  from  beta  decay.  This  rough  calculation  benchmarks  approximately 
how  many  neutrals  are  required  to  traverse  across  magnetic  field  lines  for  subsequent 
beta  decay  and  injecting  the  high-energy  electrons  into  higher  L-shells.  There  are  many 
subtleties  not  addressed  in  this  rough  calculation  such  as  the  time  dependency  of  beta 
decay  or  the  random  pitch  angle  of  high-energy  electrons  bom  of  beta  decay.  However, 
stipulating  that  approximately  15%  of  the  high-energy  electrons  born  of  beta  decay  were 
trapped  following  the  STARFISH  test  does  help  judge  that  these  results  provide  a 
reasonable  contribution  to  the  physics  underlying  the  observations  of  STARFISH. 


W.N.  Hess  “The  Atrificial  Radiation  Belt  Made  on  July  9,  1962”  The  Journal  of  Geophysical 
Research,  68  (1963)  680. 
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D,  ANALYSIS  OF  DATA 


1.  Atmospheric  profiles 

The  different  atmospherie  profiles  in  the  simulations  affeeted  the  results  of  the 
amount  of  neutrals  that  slipped  ahead  of  the  baekground  plasma.  The  eonstant 
atmosphere  resulted  in  roughly  an  order  of  magnitude  fewer  neutrals  in  front  of  the 
baekground  plasma  than  the  exponentially  decreasing  atmosphere.  This  difference  can  be 
first  explained  in  terms  of  “freeze  out.”  In  the  exponentially  decreasing  atmosphere,  the 
chance  of  either  type  of  charge  exchange  decreasing  exponentially  with  radial  distance  so 
after  the  initial  round  of  charge  exchange,  the  charge  state  remains  constant  (or  “frozen”) 
for  the  remainder  of  the  simulation.  In  a  constant  atmosphere,  the  density  of  the 
background  oxygen,  which  is  the  catalyst  for  charge  exchange,  never  changes  so  particles 
undergo  multiple  charge  exchanges  throughout  the  whole  simulation.  The  charge  state 
continues  to  alternate  and  those  debris  ions  that  become  neutral  undergo  electron  loss  and 
become  ionized  again. 

The  second  contribution  the  atmosphere  makes  to  the  effect  of  charge  exchange  is 
through  the  collisionless  shock.  A  comparison  of  Figures  23  and  24  show  that 
collisionless  shock  impacts  the  number  of  neutrals  formed.  The  transfer  of  energy  from 
the  debris  ions  to  the  background  ions  has  a  two-fold  effect  on  number  of  charge 
exchange  interactions.  First,  the  speed  of  ions  undergoing  collisionless  shock  is  decreased 
and  the  speed  directly  impacts  the  calculations  of  the  number  of  charge  exchange 
interactions  in  a  time  step.  Second,  the  cross  sections  for  both  electron  capture  and 
electron  loss  have  energy  dependence  in  them.  Collisionless  shock,  which  transfers 
energy  from  the  debris  ions  to  the  background  ions,  decreases  the  energy  used  in  the 
calculation  of  the  cross  section.  The  exponential  atmosphere  generates  less  collisionless 
shock  as  the  number  density  of  the  background  plasma  decreases  with  radius  and  this 
collisionless  shock  impacts  the  amount  of  charge  exchange. 

2,  Elemental  Considerations 

The  differences  in  the  percentages  of  neutrals  produced  between  the  various 
elements  have  several  explanations.  First,  the  atomic  differences  among  the  elements 
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results  in  different  cross  sections  (for  both  electron  capture  and  electron  loss).  The 
features  of  the  different  atomic  configurations  results  in  more  or  less  electromagnetic 
shielding  of  the  valance  electrons  affecting  the  cross  sections  for  both  electron  capture 
and  electron  loss.  Second,  differences  in  the  percentage  of  neutral  uranium  between  the 
various  runs  requires  another  explanation.  Between  the  various  simulations,  the  second 
debris  element  (lead,  zirconium,  xenon)  had  varying  number  densities  to  keep  the  overall 
energy  in  the  simulations  constant.  Varying  the  number  densities  varied  the  coupling 
coefficient  and  affected  the  path  of  the  uranium  ions.  As  the  neutral  oxygen  density  is 
related  to  the  path,  the  resulting  collisionless  shock  affected  the  percentage  of  neutral 
uranium  atoms  formed  across  the  various  simulations. 

E,  RECOMMENDATIONS  FOR  FUTURE  WORK 

The  purpose  of  this  work  was  to  determine  the  impact  of  collisionless  shock 
formation  driven  by  multiple  species  of  ions  and  charge  exchange  processes  on  the 
spatial  distribution  of  debris  from  a  HANE.  This  work  used  PIC  simulations  to  identify 
the  effect  of  charge  exchange  processes  on  previously  fully  coupled  debris  plasma 
expansions  and  showed  it  produced  significant  results  when  included  in  the  simulation. 
This  physics  should  be  included  in  a  higher  fidelity  model  (like  LLNL’s  KIM3D)  since  it 
appears  to  have  measurable  and  realistic  effects  on  the  simulation  in  simpler  geometry. 
Buffering  limitations  inherent  within  ZMR  prevented  running  more  than  two  species  at  a 
time  with  charge  exchange  enabled,  as  conserving  energy  and  momentum  resulted  in  an 
exponential  growth  in  simulation  particles  over  time.  Employing  a  different  buffering 
scheme  would  alleviate  the  geometric  growth  in  simulation  particles  that  hindered  ZMR. 
Realistically,  the  debris  from  a  HANE  initially  consists  of  thousands  of  different 
isotopes;  while  most  are  very  short  lived.  After  the  debris  has  moved  out  several 
kilometers,  it  will  still  consist  of  substantially  more  than  just  two  elements.  The  elements 
used  in  this  work  are  from  the  peaks  of  the  double  humped  curve  of  fission  fragments 
(Eigure  3)  and  also  match  the  elements  used  by  LENT  for  other  simulations.  The  use  of 
more  than  two  debris  species  will  impact  the  final  results.  As  the  shock  formation  is  non¬ 
linear  and  does  impact  the  charge  exchange,  the  inter-species  influence  will  have  impact 
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the  shock  formation.  The  differences  in  the  shock  formation  will  impact  the  charge 
exchange. 

This  work  assumed  all  the  debris  started  at  charge  state  +1.  In  the  simulation  with 
charge  exchange,  ions  could  either  be  +1  or  neutral  and  switched  back  and  forth  based  on 
both  the  processes  of  electron  capture  and  electron  loss.  Charge  states  higher  than  +1  are 
possible  and  ZMR  was  modified  to  handle  higher  charge  states.  A  higher  charge  state  ion 
would  have  a  smaller  Larmor  radius  and  equations  14  and  30  suggest  they  would  drive 
the  collisionless  shock  differently.  Starting  an  ion  at  charge  state  +2  presented  the  same 
geometric  growth  of  simulation  particles  that  plagued  multispecies  runs.  Again,  a  better 
buffering  system  would  be  able  to  handle  higher  charge  states.  While  ZMR  can  handle 
higher  charge  states,  there  is  no  current  mechanism  to  handle  charge  exchange  in  integer 
values  greater  than  one;  that  is,  the  code  can  not  yet  handle  multiple  electron  losses  or 
captures  in  a  single  event.  The  current  limit  remains  changes  to  the  charge  state  of  ±1 
only. 

More  robust  codes  could  better  account  for  the  real  atmosphere.  This  work  used 
both  a  constant  density  atmosphere  and  an  exponentially  decreasing  atmosphere.  The 
angular  symmetry  of  ZMR  makes  neither  of  these  assumptions  physical.  The  real 
atmosphere  is  exponentially  decreasing  in  the  vertical  direction  and  exponentially 
increasing  in  the  downward  direction.  It  is  relatively  constant  parallel  to  the  Earth.  From 
the  curvature  of  the  L-shells,  a  bound  particle  will  encounter  a  variety  of  densities  as  it 
travels.  A  charged  particle  bound  to  a  field  line  of  the  magnetic  field  will  encounter  a 
different  atmospheric  density  at  the  magnetic  equator  and  at  its  mirror  points.  None  of 
those  conditions  were  included  in  this  simulation  and  should  be  examined  in  future  work 
to  more  accurate  replicate  the  physical  geometry  as  the  ultimate  goal  is  to  predict  the 
locations  and  charge  state  of  the  debris  from  any  HANE. 

The  importance  of  the  neutral  fraction  Y  (equation  49)  bears  further  examination, 
since  this  quantity  impacts  the  total  number  of  neutral  particles.  It  suggests  those  neutrals 
with  high  ionization  potentials  (like  the  Group  7  and  8  elements  from  the  periodic  table) 
would  have  the  largest  impact  on  the  total  number  of  neutrals.  Further  code  development 
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that  focused  on  those  elements  that  had  the  highest  probability  of  forming  neutrals  would 
be  an  effieient  method  to  draw  out  the  physieal  effeets  of  those  elements. 

E,  CONCLUSIONS 

This  dissertation  examined  the  motion  of  the  debris  from  a  high  altitude  nuelear 
explosion  (HANE)  like  STARFISH  to  determine  if  the  physies  of  eollisionless  shook 
generated  by  multiple  speoies  and  oharge  exohange  supported  the  observation  of  high- 
energy  eleotrons  at  higher  L-shells.  It  was  shown  that  the  inolusion  of  oharge  exohange 
does  impaot  the  kinetio  motion  of  the  debris  ions  and  can  move  them  to  looations  not 
previously  predioted  but  aligned  with  the  measured  data.  Charge  exohange  aocounts  for 
previously  unreported  movement,  whioh  positions  portions  of  the  debris  to  looations 
where  later  beta  deoay  would  injeot  high-energy  eleotrons  oonsistent  with  observations. 
The  effeets  of  multiple  species  in  the  debris  plasma  is  also  found  to  be  signifioant  and 
demonstrates  the  importanoe  of  the  non-linearity  of  the  shook  formation  on  the  oharge 
state  of  the  debris  ions  while  employing  the  charge  exohange  prooesses.  Inolusion  of  the 
eleotron  oapture  and  eleotron  loss  oomputational  routines  in  a  higher  fidelity  oode  is 
important  to  get  relevant  physieal  results.  If  implemented,  obtaining  the  proper  results  is 
likely  to  be  well  worth  the  additional  oomputational  time  beoause  this  work  has  shown 
the  physieal  results  to  be  signifioantly  impaoted  by  the  inolusion  of  oharge  exohange. 
Collisionless  shook  driven  by  multiple  speoies  was  also  modeled  and  did  not  prove 
significant  as  a  stand-alone  hypothesis,  but  the  importance  of  eollisionless  shook  on 
oharge  exchange  was  demonstrated. 

An  overview  of  the  simulations  examining  oharge  exohange  is  illustrated  in 
Figure  41.  After  the  nuelear  detonation  (red  star  in  Figure  41),  the  debris  streams 
perpendicular  to  the  ambient  magnetio  field.  The  interaotion  of  this  debris  with  the 
baokground  plasma  of  flash  ionized  monatomio  oxygen,  shown  as  two  dashed  blue  lines 
representing  both  density  profiles,  generates  the  eollisionless  shook.  As  the  debris 
oontinues,  it  internets  with  the  neutral  monatomio  oxygen  in  the  baokground  and 
undergoes  oharge  exohange.  Onoe  neutral,  the  partioles  oontinue  on  their  trajeotory 
unaffeoted  by  the  eleotromagnetio  foroes.  After  several  seoonds  of  travel,  those  neutrals 
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were  found  to  arrive  at  much  higher  L-shells  and  undergo  beta  decay  introducing  high- 
energy  electrons  similar  to  those  observed  in  STARFISH. 

Ambient 

magnetic 

field 


★ 


injecting  high-energy 
electrons  to  radiation  belt 


No  charge 
exchange  within 
first  50  km 


Exponential  0*°  profile  as 
catalyst  for  charge  exchange 


Figure  41.  Graphic  summary  of  simulations  examining  charge  exchange. 

Charge  exchange  is  shown  to  be  a  significant  contributing  factor  to  this  process 
and  explains  the  observation  of  high-energy  electrons  at  higher  L-shells  following  a 
HANE;  the  physics  of  charge  exchange  should  be  included  in  future  high  fidelity 
modeling.  While  collisionless  shock  did  not  have  an  independent  impact  on  the  motion  of 
the  debris,  when  included  in  conjunction  with  charge  exchange  it  did  affect  the 
simulation  outcomes  and  future  simulations  should  employ  both  principles  to  capture  the 
relevant  physics. 
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In  the  course  of  this  research,  both  electron  capture  and  electron  loss 
computational  subroutines  were  developed  and  implemented  in  a  particle-in-cell  hybrid 
code.  The  quantifiable  differences  of  including  these  subroutines  was  determined  and 
proved  to  be  important. 

While  this  work  was  not  intended  to  replicated  the  precise  results  from 
STARFISH,  it  did  prove  the  value  of  including  the  physics  of  collisionless  shock  and 
charge  exchange  into  simulations  in  order  to  examine  the  motion  of  the  debris  across  the 
magnetic  field  lines  to  higher  L-shells  where  they  beta  decay  and  inject  high-energy 
electrons. 
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